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Abstract 



This thesis consist of two parts, post-Newtonian modes of relativistic sys- 
tems and the stabihty curve of r-modes of neutron stars. In part I, we use 
the post-Newtonian (pn) order of Liouville's equation to study the normal 
modes of oscillation of a spherically symmetric relativistic system. Pertur- 
bations that are neutral in Newtonian approximation develop into a new 
sequence of normal modes. In the first pn order; a) their frequency is an 
order q smaller than the classical frequencies, where q is a. pn expansion 
parameter; b) they are not damped, for there is no gravitational wave ra- 
diation in this order; c) they are not coupled with the classical modes in q 
order; d) because of the spherical symmetry of the underlying equilibrium 
configuration they are designated by a pair of angular momentum eigen- 
numbers, {j,m), of a pair of phase space angular momentum operators 
{j'^,Jz)- The eigenmodes are, however, m-independent. Hydrodynamics of 
these modes is also investigated; a) they generate oscillating macroscopic 
toroidal motions that are neutral in classical case ; and b) they give rise to 
an oscillatory g^oi component of the metric tensor that otherwise is zero in 
the unperturbed system. The conventional classical modes, which in their 
hydrodynamic behavior emerge as p and g modes are, of course, perturbed 
to order q. These, however, have not been of concern in this work. 

In part 11, stability curve of r-modes of neutron stars are calculated by 
considering vorticity-shear viscosity coupling. The coupling is predicted by 
kinetic theory, a causal theory of fluid rather than the Navier-Stocks theory. 
We calculate this coupling and show that it can in principle significantly 
modify the stability diagram at lower temperatures. As a result, colder 
stars can remain stable at higher spin rates. 



As an application, the loss of angular momentum through gravitational 
radiation, driven by the excitation of r-modes, is considered in neutron stars 
having rotation frequencies smaller than the associated critical frequency. 
We find that for reasonable values of the initial amplitudes of such pulsation 
modes of the star, being excited at the event of a glitch in a pulsar, the total 
post-ghtch losses correspond to a negligible fraction of the initial rise of the 
spin frequency in the case of Vela and older pulsars. However, for the Crab 
pulsar the same effect would result, within a few months, in a decrease in 
its spin frequency by an amount larger than its glitch-induced frequency in- 
crease. This could provide an explanation for the peculiar behavior observed 
in the post-glitch relaxations of the Crab. 
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Part I 



Post-Newtonian modes 
relativistic systems 



Chapter 1 
Introduction 



Chandrasekhar's formulation of post-Newtonian (pn) hydrodynamics is 
among the pioneering ones. He generahzed Eulerian equations of Newtonian 
hydrodynamics to pn order consistent with Einstein's field equations, and 
applied them to obtain the pn corrections to the equilibrium and stability 
of uniformly rotating homogeneous masses. Blanchet, Damour and Schafer 
[^] studied the gravitational wave generation of a self gravitating fluid by 
adding an appropriate term to pn equation of hydrodynamics. Cutler 
employed the pn hydrodynamics and a perturbation technique to derive 
an expression for the pn correction to Newtonian eigenfrequencies. Cutler 
and Lindblom adopted Cutler's method to calculate numerically the 
oscillation frequencies of the / = m /-modes of rapidly rotating polytropic 
neutron stars. 

In this work we study normal modes of a non-rotating relativistic system 
in pn approximation through the relativistic Liouville's equation rather than 
the relativistic hydrodynamics. The reason for doing so is to avoid ther- 
modynamic concepts being incorporated into hydrodynamics. Liouville's 
equation is a purely dynamical theory and free from such complexities. 
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Furthermore in many cosmological and astrophysical situations, an ideal- 
ized fluid model of matter is inappropriate, and a self-consistent microscopic 
model based on relativistic kinetic theory ^ gives a more detailed physical 
description. A well-known example is the case of coUisionless particles, as 



for cosmological neutrinos and photons ||T^, or stellar clusters in equilib- 
rT|. Also there are other examples among non-equilibrium evolving 



rmm 



systems, such as stellar clusters with collisions |jT2[, the early evolution of 
a FRW universe into an anisotropic Bianchi universe or an inhomogeneous 
universe via a disturbance of the equilibrium collisional balance The 

I, or mixture of 



relativistic transport of photons and cosmic rays 
cosmic elementary particles ||T6[ are other non-equilibrium situations suited 
to a kinetic approach. 

The kinetic theory offers an alternative approach to describe the matter, 
rather than the phenomenological fluid dynamics and its associated ther- 
modynamics. For example, the standard thermodynamics of fluids violates 
causality and is unstable A casual and stable generalization emerges 



from the kinetic theory, as developed by Israel and Stewart |T8|. In some 
cases, a fluid model leads to the loss of information and is unable to account 
for certain effects. For example. Landau-type damping of gravitational per- 



turbations by a kinetic gas is not present in the fluid models [0]. Another 
example is the rotational perturbations coexisting with initial singularity, 
which is impossible in the fluid case pO[ . 



In compiling this work we have relied heavily on the following studies 
dealing with various aspects of Liouville's, Liouville-Poisson's and Antonov's 
equations. 

0(3) symmetry and mode classiflcation of classical Liouville's equation 



for spherically symmetric potentials was studied by Sobouti [|2T|. Simple 
harmonic potentials in one, two, and three dimensions were discussed by 
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Sobouti |2^. He obtained exact and complete eigensolutions by means of 
raising and lowering ladders for Liouville operator. Furthermore, he investi- 
gated potentials of self gravitating spheres, oblate or prolate spheroids, and 
ellipsoids in details. A systematic method to elaborate the symmetries of 
Liouville's equation for an arbitrary potential were introduced by Sobouti 
and Dehghani ||2^. They showed that the symmetry group of potential 
is GL(3, c) and classified eigenmodes of Liouville's equation for quadratic 
potential. 0(4) symmetry of potential were obtained by Dehgahni and 
Sobouti ||2^. Dynamical symmetry of Liouville's equation for potential 
was worked out by Dehghani and Sobouti Dynamical symmetry group 
of general relativistic Liouville's equation was discussed by Dehghani and 
Rezania p6[. In particular they found that in de Sitter's space-time the 
group is S0(4,l) (g) S0(4,l). 

In applications to self-gravitating systems the pioneering work was done 
by Antonov |^ . He reduced the linearized Liouville- Poisson equations to a 
self adjoint operator in phase space. Further elaborations on Antonov's 
equation were made by Lynden-Bell ||28|, Milder pOf, Lynden-Bell and 
Sanitt Ipser and Thorne These authors were concerned with the 
stability of a given isotropic distribution function. Stabilities of anisotropic 
distribution functions were investigated by Doremus et. al. Dore- 
mus and Feix |3^, Gillon et. al. [^, Kandrup and Sygnet |Q. Attempts 
to solve the linearized Liouville-Poisson equation for eigenfrequencies and 



eigenmodes of oscillations were made by Sobouti [pTj pq , |39| . Further and 
more transparent exposition of mode classification and mode calculations 
were given by Sobouti and Samimi EO] and Samimi and Sobouti El . 



Here, using the standard pn expansion of the metric components |42 



we derive the pn approximation of Liouville's equation {pnl). In the time- 
independent case, we show that a generalization of classical integrals (en- 
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ergy and angular momentum, say) are the static solutions of pnl. In time- 
dependent regimes, the effect of the pn corrections on the known solutions 
of the classical equation can be analyzed by the usual perturbation tech- 
niques. Whatever the procedure, the first order corrections on the known 
modes will be small and will not change their nature. We will not pursue 
such issues here. The main interest of this work is to study a new class of 
solutions of pnl that originate solely from the pn terms and have no prece- 
dence in classical theories. It is not difficult to anticipate the existence of 
such modes. Perturbations on an equilibrium state, that are functions of 
classical integrals do not disturb the equilibrium of the system at classical 
level. That is they do not induce restoring forces in the system. They, 
however, do so in the pn regime, and make the system oscillate about the 
pn equilibrium state. Such perturbations may be considered as a class of 
infinitely degenerate zero frequency modes of the classical system. The pn 
forces unfold this degeneracy and turn them into a sequence of non zero fre- 
quency modes distinct and uncoupled from the other classical modes. We 
have termed them as pn modes. 

A hydrodynamic analog of pn modes is the following. In spherically 
symmetric fluids, toroidal motions are neutral. Sliding one spherical shell 
over the other is not opposed by a restoring force. A small magnetic field 
or a slow rotation (mainly through Coriolis forces) gives rigidity to the 
system. The fluid resists against such displacements and a sequence of 
well defined toroidal modes of oscillation develop. See Sobouti Hasan 
and Sobouti ||^, Nasiri and Sobouti and Nasiri for examples and 
typical calculations. 

The plan of this part is as follow. In chapter 2 we briefly review Liou- 
ville's equation both in classical and relativistic regimes. We introduce a 
distribution function and its equation of evolution. Macroscopic quantities 
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associated to the distribution function are discussed. 

In chapter 3 we adopt the post-Newtonian (pn) approximation to study a 
self gravitating system imbedded in an otherwise flat space-time. We obtain 
the pn approximation of Liouville equation (jml). We find two integrals of 
pnl that are the pn generalizations of the energy and angular momentum 
integrals of the classical Liouville's equation. Post-Newtonian polytropes, 
as simultaneous solutions of pnl and Einstein's equations, are discussed and 
calculated 

In chapter 4 we give the pn order of the linearized Liouville equation 
that governs the evolution of small perturbations from an equilibrium state. 
We extract the equation for a sequence of new modes that are generated 
solely by pn force but are absent in classical regime. We explore the 0(3) 
symmetry of the modes and classify them on basis of this symmetry. We 
study hydrodynamics of these modes. We seek a variational approach to 
the calculation of pn modes and give numerical values for polytropes. 

Post- Newtonian approximation is reviwed in appendix A. In appendix 
B coordinates transformation that we need to extract Liouville's equation 
in pn approximation, is discussed. In appendix C post-Newtonian hydro- 
dynamics arc recovered by integration of pnl over v-spacc. Simultaneous 
eigensolutions of and operators are constructed and elaborated in ap- 
pendix D. 
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Chapter 2 

Liouville's equation 
Classical and Relativistic 



Kinetic theory has expanded in classical, quantum, and relativistic direc- 
tions [0. Classical kinetic theory is the foundation of fluid dynamics and 



thus is important to aerospace, mechanical, and chemical engineering. It 
is also relevant to many problems in astrophysics, for example the stabil- 
ity and evolution of stellar systems. Quantum kinetic theory is applicable 
to problems in particles transportation, radiation through material media, 
etc., which are important in solid state and laser physics. Relativistic ki- 
netic theory has became important in certain plasma physics. It is also 
used to study the evolution of relativistic stellar systems and the dynamics 
of cosmological fluids. 

In its most elementary version the kinetic theory of a simple gas re- 
lies on the concept of N pointlike particles which may interact with each 
other. Collisions are assumed to establish a local or global equilibrium of 
the system. Between the collisions the particles move on geodesies of a given 
spacetime. 
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Technically, the gas particles are described by an invariant one-particle 
distribution function governed by Boltzmann or Liouville's equation. The 
latter is best applicable to dilute gases where collisions may be neglected. 
The macroscopic fluid dynamics for such a system may be obtained in terms 
of the first and second moments of the distribution function. A gas, how- 
ever, is the only system for which the correspondence between microscopic 
variables, governed by a distribution function, and phenomenological fluid 
quantities is sufficiently well understood. 

The goal of this chapter is to introduce the Liouville's equation both 
in classical and relativistic cases. Liouville's equation gives the time evo- 
lution of probability distribution function. It provides the dynamical basis 
of statistical mechanics , both at and away from equilibrium . Its solution 
enables one to calculate the ensemble average of any dynamical quantity . 

In section 2.1 classical Liouville's equation is discussed. For the apph- 
cation in astrophysical problems, the linearized Liouville-Poisson equation 
is introduced. Relativistic Liouville's equation is considered in section 2.2. 
In section 2.3, macroscopic quantities are discussed. 

2.1 Classical Liouville's equation 

For a system of degrees of freedom, phase space is a 2A^ dimensional 
space whose axes are the {xi,Pi) variables. Thus the state of the system 
at any given instant (x, p) is a single point, which is usually called system 
point, in 2N dimensional phase space. Time is exhibited explicitly in 2A^-|-1 
dimensional phase space, a 2N phase space with an additional orthogonal 
time axis. As time evolves, the system point moves on a system trajectory, 
[x(t),p(t)], which is a curve in 2A^ + 1 space. The system trajectory is 
determined by solving the equations of motion: 
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dp 

P = ~. (21b) 
ax 

where H is the Hamiltonian of the system. It is clear that the state of the 
system will be specified uniquely by 2N initial constants, (x(0),p(0)). 

An abstract collection of a large number of independently identical sys- 
tem points is called an ensemble. An important property of the ensemble is 
that trajectories of the ensemble can never cross in phase space. This follows 
from the fact that for a system with N degrees of freedom the system of 
trajectory, Eq. (|2.1| ), is uniqely specified by 2N initial values, [x(0),p(0)]. 

Consider an infinitesimal volume in phase space surrounding a given 
system point at time t = 0. In the course of time the system points defining 
a volume element move about in phase space and the volume contained by 
them will take on different shape as time progresses. The Liouville theorem, 
however, states that the size of a volume element in the phase space remains 
constant under canonical transformations induced by the Hamiltonian, i.e. 
the Jacobian of a canonical transformation is unity |^ . 

Let dAf denotes the number of system points in a phase space volume 
element. It remains constant. For, a system point initially inside can never 
get out, and one outside can never enter the volume. Indeed, if some system 
point were to cross the border, its trajectory would intersect a trajectory 
of a system point defining the boundary surface. But this is not possible. 
For, if two trajectories were to coincide at one time, they would coincide at 
all the times. Hence the number of the system points, d-JV, within a volume 
element of phase space, dT, remains constant. In other words the probability 
density, /(x, p,t) = J\f^^dJ\f/dT, should remain constant in time. That is 



9 



1 = 0. (2.2a) 
This is the Liouville's equation. Assuming / is differentiable, we obtain 

df df df . df . ^ 

J- = -jL _| ■Lj.. _| ■Lp. = 0. (2.2b) 

dt dt dxi dpi 

Taking equations of motion, Eq. (pTip, into account, Eq. ( |2.2b| ) becomes 

dl^dl^mdl_mdl^^ (2 3) 
dt dt dpi dxi dxi dpi 

It is convenient to write it as 



^|^ = £^7, (2.4a) 
where classical Liouville's operator, is the linear operator 

_ -^dH d OH d 

dpi dxi dxi dpi 

As it will be shown later, the reason for including i is to render C^^ Hermi- 
tian. In terms of the formal solution of Eqs. (p.4|) is 



/(x,p,t)=e-*^V(x,p,0). (2.5) 

It is easy to show that if the initial /(x, p,0) is an acceptable distribution 
function, /(x, p, t) will be an acceptable one at all later time. In particular 

/(x,p,t)>0; Vt, (2.6a) 
j f{x,p,t)dT = 1; Vt. (2.6b) 



See Balescu 48 
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In this section, we obtained classical Liouville's equation in general form. 
To solve the equation for specific problem, we must first define the Hamil- 
tonian of the system. 
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2.1.1 Properties of C 

In this section we review some important properties of C^K For more details 
see [0, |2g, H, 11, |T], H . 



The Hilbert space: An axiomatic study of the eigensolutions of classical 
Liouville's equation requires introduction of a Hilbert space. A Hilbert 
space, H, is defined to be the space of complex square integrable functions 
of phase coordinates (x, p) that vanish at the phase space boundary of the 
system: 

H : /(x, p); J f*fdT = finite, /(boundary) = . (2.7) 

Integrations in Ti. are over the volume of the phase space available to the 
system. 

Hermiticity: C^^ is Hermitian in TC, i.e. 

J g*{C^^f) dV = j {C^'gYf dV; gJeH. (2.8) 

This is proved by integrating by ( p.8| ) by parts and letting the integrated 
terms vanish at boundary. 

Real eigenfrequency: The eigenfunctions /n(x, p) and eigenvalues Un are 
defined by 

C''fn = UJnfn. (2.9) 

Hermiticity of ensures that are real and eigenfunctions belonging 
to distinct eigenvalues are orthogonal. 
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Completeness of eigensolutions: We can further impose the normahza- 
tion condition on /„'s, 

I fUndT = 6mn, (2.10) 

and obtain an orthonormal set. We shall also assume that they also form a 
complete set. 

Since classical Liouville operator is purely imaginary, C^^* = —C^'', its eigen- 
solutions have following properties: 

(a) Eigensolutions belonging to non-zero eigenvalues are complex, ie. 

/(x,p) = m(x,p) +if(x,p) . (2.11) 

(b) If {uj, f) is an eigensolution, {—uj, /*) is another eigensolution. 

(c) /*/ is an integral of motion, ie. C'^\f*f) = 0. 

(d) [{n — 171)00, f*"^f"'] is an eigensolution with n,m = positive integer. 

(e) Eigenfunctions belonging to non-zero eigenvalues integrate to zero: 

JfdT = 0, uj^O. (2.12) 
2.1.2 Linearized Liouville-Poisson's equation 



In applications to astrophysical problems, many investigators |2^-[^, have 
often used the linearized Liouville-Poisson equation to study stability of the 
perturbed a self-gravitating stellar system. 



In this section, we follow Sobouti 37, 38, 3£], Sobouti and Samimi 40 



and Samimi and Sobouti |^ to introduce the classical linearized Liouville 
equation. 
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For a coUisionless self-gravitating stellar system the classical Liouville's 
equation, Eqs. ( |2.4|) , for distribution function, F(x, p,t), becomes 



t— = C^^F, (2.13a) 
rci ■( d dU d\ 



' dxi dxi dpi J 

where the potential f/(x, t) is the solution of Poisson's equation 

[/(x, t) = -G J F(x , p', t) I X - X dT' . (2.14) 

The Hamiltonain used here is the energy of the system, E = + f/(x, t). 
It is easy to see that the energy is an integral of C'^'' in an equilibrium state. 
Furthermore, for spherically symmetric potentials, the angular momentum, 
h = ^ijkXjPk, is another integral of C^K 

To find the linearized equation, let F F{E) +(5F(x, p, t), where F(E) 
is an equilibrium distribution function, and | 6F \ < F{E) for all (x, p, t) 
is a perturbation on F{E). Accordingly, the potential splits into a large and 
small terms, U{r) + 6U{'x,t) where r =| x |. Substituting in Eqs. ( |2.13aD 
and ( p.l4|) , in the first order we find 



. dSF _ ^^^^ ^ .dFdSU 
dt dpi dxi 

= C'^SF + GFeC"^ J 5F(x', p', t) I X - x' dV , (2.15) 

5f/(x, t) = -G j SF{x.', p', t) I X - x' dT' , (2.16) 

where is now constructed with the time-independent potential f/(r) and 
Fe = dF/dE. 
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Let 5F =\Fe 1^/2 /(x, p, t), see Sobouti g^, then Eq. (glsD can be written 
as 

■^f Af 

(2.17) 

^/ = /:^7+ Gsign(Fs)|Fs|i/2/:^'y' |F^|V2/(x',p',t)|x-x'ritir'. 

It is easy to show that for the hnearized equation the classical energy, E = 
+ U, is not an integral, but the angular momentum is. Conservation of 
angular momentum means that the operator A has 0(3) symmetry. 
Let / = /_(x, p,t) + i/+(x, p,t), where /_ and /+ are odd and even in p, 
respectively. Substituting this in Eq. ( p.l7|) , and decomposing it into odd 
and even components, we find 



AU, (2.18a) 



dt 



Af- , (2.18b) 



dt 

where A is odd in p. Eliminating /+ we obtain a wave equation for /_ 



dt 
where 



=^V-, (2.19) 



A^f^ = /:^'V-+G'sign(F^) I Fe C^'l I Fe \'^' /:^'7-(x', p', t) | x-x' \-' 

(2.19a) 

Equations ( |2.19| ) are Antonov's equation. /_ and /+, calculated from Eqs. 
(|2.19|) and (|2.18b| ), give a solution of the linearized Liouvolle-Poisson equa- 



tion. Assuming sinusoidal time dependence for /_(x, p,t) = /_(x, p)e*'^*, 
we find 
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A'f. = uj'U ■ (2.20) 

Equation (|2.2CI|) is an eigenvalue problem with eigenfrequancy u. Sobouti 
and Samimi ||4^ proved that the A is not Hermitian in 7i, i.e. A ^ . 
However, by decomposing the Hilbert space H into odd and even subspaces 
in p, they showed that A^ is Hermitian on odd subspace: 



J flA'UdT = J I C'^'f. \UT+ G sign{FE) x 

X J \ Fe C^'U I C'^'f'_ I X - x' |-i dTdV = real. 

(2.21) 



Equation ( |2.21| ) ensures that the eigenfrequancies are real. 

0(3) symmetry of A: For spherically symmetric potentials, the invariance 

of A under rotation of both x and p coordinates was established by Sobouti 



and Samimi |^0[. The corresponding angular momentum operator in phase 
space is 



with the angular momentum algebra 



■' h u-' — 



(2.22) 







(2.23a) 
(2.23b) 



We note that Jj rotates simultaneously both x and p coordinates. Commu- 
tation of Ji with C^^ was first proved by Sobouti pT| : 
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(2.24) 



Sobouti and Samimi |^ extended the same to A, 



[A, Ji] = (2.25) 

An important consequence of Eqs. (|2.23| ) and ( p.25|) is the mutual commu- 
tation of the following set of operators 

J^ jj = 0. (2.26) 

The implication of Eq. ( p. 261) is clear. The eigenfunctions of A"^ can simul- 
taneously be the eigenfunctions of and J^. In other words, the eigen- 
functions of A'^ can be classified into classes specified by eigennumbers j, 
m of and J?. See for more details [EO . 



2.2 Relativistic Liouville's equation 

In section 2.1, we introduced classical Liouville and the linearized Liouville- 
Poisson equations. We reviewed some properties of these equations, that 
help one to extract their eigenfunctions and eigenfrequencies. The goal of 
the present section is to introduce the distribution function and its equation 
of evolution in general relativity. For this purpose we need to introduce the 
pahse space on which such a function is defined. 

2.2.1 Distribution function 

Consider a single test particle with mass m which moves in a gravitationally 
curve spacetime. Its motion is determined by the geodesic equation 

^ rfA' dX - dX ^ '-'^^ ^ ^' ^ ^' 

where A is an afiine parameter defined by the requirement that p'^ be the 
4-momentum. Hereafter, F^^ are Christoffel symbols associated with the 
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metric Qfj^^. Note that if there are non gravitational forces (e.g. electromag- 
netic forces) then we have to modify this equation. 
The rest mass of the particle is defined as 



= -p^p^. (2.28) 



Thus, according to Eq. ( p. 271) , the state of the particle is determined by the 



pair {x^,p^). The phase space is then the tangent bundle over the spacetime 
manifold, i.e. 

T = {{xf',pf'),x^' e M^p^" (2.29) 

where Ai is the space-time and 7^ is the tangent space to Ai at x'^. From 
now on, Greek indices run from to 3 and Latin indices run from 1 to 3. 

The volume element on % supported by the displacements dpi, dp2, dp^, dp^ 
(with components dj?" etc.) is 

Hp") = e^p^sdpyp^2dpldpi (2.30) 

where ex^/s^ is the totally antisymmetric tensor such that eoi23 = \^—g- We 
also define 7r+(p'^), the volume element corresponding to the subspace of 7^ 
such that is non-spacelike and future directed, 

iT+ip'') = H{-pX)H{-p'Mp'^), (2.31) 

where H is the heavyside step function 



H{x) 



1 if X > 0, 
otherwise. 



and an arbitrary timelike vector field. 

Tx is sliced in hypersurfaces, Vm, of constant m called the mass-shell, and 
defined by 

VU^n = {p" e r,,p% = -m^^p^u^ > O} . (2.32) 
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The volume element of Eq. (p. 301) on T can then be decomposed into a 
volume element, rmTm, on Vm by 

7r+(j9^) = m7r,n{p^)dm. (2.33) 

The factor m allows one to include particles of zero rest mass (see Ehlers 0]). 
This defines the induced volume element mTTrn{p^) on Vm- If we introduce 
an arbitrary future directed unit timelike vector (i.e. satisfying u^u^ = 
— 1), the 3- volume supported by the three displacements dxi, dx2, dx^ (with 
components dx"^ etc.) in the hypersurface perpendicular to u'^ is 

dV{u^) = exaf^^u^dxyx^dxl (2.34) 

We now consider a single fluid composed of particles of all masses. The 
distribution function, f{x^^p^) will be defined as the mean number of par- 
ticles (on a statistical set) in a volume dV around x^ and 7r(]9^) around 
measured by an observer with 4-velocity u^, 

dN{xf',pf') = /(a;^,p^)(-p%)rfy(M'^)7r(p'^). (2.35) 

The assumptions involved in its existence have been discussed in details by 
Ehlers [Q. Synge has demonstrated that {—p^u^)dV{u^) is indepen- 



dent of M^. This implies that the distribution function is a scalar. Moreover, 
/(x^,p^) > for all x'^ and all allowed p^. 

For a gas, dN is the number of particles in a volume dVTT{p^) thus the 
smoothness of / depends on the existence of a sufficient number of particles. 

2.2.2 Relativistic Liouville opertaor 

The equations of motion, Eq. ( |2.27|) , define on T the Liouville operator, 

„ d dp'^ d d , , 

^x^' dX dpi" dX ^ ^ 



which characterizes the rate of change of / along the particle's worldlines. 
Using (p.27|), this operator can be rewritten as 



d 



rt>v 



d 



(2.37) 



The fact that the mass m of Eq. ( p.28|) is a scalar constant on each phase 
orbit leads to 

L{m^) = 0. (2.38) 

The Boltzmann equation states that this rate of change is equal to the rate 
of change due to collisions, i.e. that 



Cf = C[f]. 



(2.39) 



C[f] is the collision term and encodes the information about the interac- 
tions between the particles of the fluid. 



If we now consider a system of fluids (labelled by each of which 

is described by its distribution function fi{x'^,p'^), the Boltzmann equation 
for a given fluid i becomes 

Cf. = T.CAf^J,]^Cm, (2-40) 
j 

Cj[fi, fj] is the collision term describing the interaction between the fluid i 
and the fluid j. For elastic collisions, it must satisfy the symmetry 

CAf^Jj] = C^[f., a (2.41) 

which means that in a collision between i and j the two distribution func- 
tions undergo the same change. If we assume the gas is dilute (i. e. the 
mean free path of particles is much greater than the range of interactions 
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between them) such that we can neglect colhsions between its particles, the 
RHS of Eq. ( p.39 ) will be zero. Therefore we find 



Cf = 0. (2.42) 

This is the general relativistic Liouville's equation. The equation describes 
the evolution of distribution function, /, of a coUisonless gas. In the mass- 
shell hypersurface, Vm, which is defined in Eq. ( |2.32| ), Liouville's equation 
for all particles with the constant mass m reduces to 

^■"^=(^a^-W|.)/ = 0. (2.43) 

where Cm is the Liouville operator, C, in the mass-shell Vm- 



2.3 Macroscopic quantities 

In section 2.2, we have introduced distribution functions and the basic evo- 
lution equations of the system in general relativity. The goal of this section 
is to define a set of macroscopic quantities from the distribution function 
and the collision term and then find the relations between these quantities. 
Given a distribution function /, at any point x'^, one can introduce, follow- 



ing Ellis et al. a set of macroscopic quantities associated with fluid 

by 



mVK..p^"Mx^,p^)7rm, (2.44) 



Vn 



where m is the mass of the particles (defined by Eq. ( |2.28| )) and a an integer. 
If the particles of a given fluid have different rest mass (this is the case e.g. 
when one is dealing with a fluid of stars or of galaxies) the above equation 
should be modified by an integration over m (see Uzan [^). We assume 
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that each distribution function vanishes at infinity on the mass shell rapidly 
enough so that all these integrals converge. 

Among all these quantities, some are important in many applications, 

n^^Xo^ = |p'^/(a;^p'^)7^^, (2.45) 
j,,u ^ ^ I (2.46) 

The vector is the number fiux vector which is used to define the average 
number fiux velocity vector and the proper density n measured by an 
observer comoving with the fiuid by 

= nv^", = -1. (2.46) 

T^jy is the energy-momentum tensor. In terms of the timelike unit vector 
field u'^, chosen as time direction, we can split the energy-momentum tensor 
under the general form 

T^,u = (p + P)Uf,Uy + Pg^,^ + 2(g^M^ + q^u^,) + vr^^, (2.47) 

the quantities p, P, g'^ and n^y being defined as 

p = T^uuV, (2.48) 

P = ^T^.h^\ (2.49) 

= -h'^n^u'', (2.50) 

vr^. = KKT^p - Ph^,, (2.51) 

where h^^i, = g^^ + M^Wy. This decomposition is the most general splitting 
with respect to the arbitrary vector field of a tensor of rank 2. The four 
quantities p, P, g^, and T^^y are respectively called the energy density, the 
pressure, the energy fiux vector, and the anisotropic stress tensor. For the 
latters one can verify, from Eqs. iWM) -( ^311) 
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qX = K = ^M-"'' = 0- (2-52) 
The energy momentum tensor, (|2.46|) , appears in the RHS of Einstein's 
field equations 



^Qiiv^ — ^4 T^iivi (2.53) 



1 o SttG, 

where R^y and R are the Ricci tensor and the scalar curvature respectively. 
Therefore, in applications to self gravitating stars and stellar systems, one 
should combine Einstein's field equations and Liouville's equation, ( |2.43| ). 
The resulting nonlinear equations can be solved in certain approximations. 

In the next chapter we adopt the post-Newtonian approximation to 
study a self gravitating system. We derive Liouville's equation in this ap- 
proximation by using the standard post-Newtonian expansion. 
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Chapter 3 

Liouville's equation in post 
Newtonian approximation 



Solutions of general relativistic Liouville's equation (grl) in a prescribed 
space-time have been considered by some investigators. Most authors have 
sought its solutions as functions of the constants of motion, generated by 
Killing vectors of the space-time in question. See for example Ehlers | [5U| , 



Ray and Zimmerman |Q, Mansouri and Rakei |^5[, Ellis, Matraverse and 



Treciokas [52], Maartens and Maharaj [Q, Maharaj and Maartens [57|, 
Maharaj |Q, and Dehghani and Rezania p6| . 

In application to self gravitating stars and stellar systems, however, one 
should combine Einstein's field equations and grl. The resulting nonlinear 
equations can be solved in certain approximations. Two such methods are 
available; the post-Newtonian (pn) approximation and the weak-field one. 
In this chapter we adopt the first approach to study a self gravitating system 



imbeded in an otherwise flat space-time. In section |3.1|, we derive the pn 



approximation of the Liouville equation (pnl). In section 3^ we find two 



integrals of pnl that are the pn generalizations of the energy and angular 
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momentum integrals of the classical Liouville's equation. Post-Newtonian 
polytropes, as simultaneous solutions of pnl and Einstein's equations, are 
discussed and calculated in section |373| . 

The main objective of this chapter, however, is to set the stage for the 
chapter 4. There, we study a class of non static oscillatory solutions of pnl, 
which in their hydrodynamical behavior are different from the conventional 
p and g modes of the system. They are a class of toroidal motions driven 
by pn force terms and are accompanied by oscillatory variations of certain 
components of the space-time metric. 

3.1 Liouville's equation in post-Newtonian 
approximation, General 

The one particle distribution function of a gas of collisionless particles with 
identical mass m, in the restricted seven dimensional phase space 

P(m): g^,U^U'' = -^ (3.1) 

satisfies grl: 

CuF = ([/-A _ r^M^U^^^)F{x^, U^) = 0, (3.2) 

where {x'^, If^) is the set of configuration and velocity coordinates in P{m), 
F{x^, f/') is a distribution function, Cu is Liouville's operator in the (x'", f/*) 
coordinates, F^^, are Christoffel's symbols, and c is the speed of light. Greek 
indices run from to 3 and Latin indices from 1 to 3. The four- velocity of 
the particle and its classical velocity are related as 

f/M = ifO^f^. = = dx'/dt), (3.3) 

where U^{x^,v^) is to be determined from Eq. (3.1). In pn approximation, 
we need an expansion of Cu up to the order (w/c)^, where w is a typical 
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Newtonian speed. To achieve this goal we transform (x^,U^) to {x^,v'^). 
Liouville's operator transforms as 

where dv^/dx^ and dv^/dlf^ are determined from the inverse of the trans- 
formation matrix (see appendix B). Thus, 



dW Q 



dxt" 2Q dx^" 
dv^ 1 



(3.5b) 



--^(t/°"' + Y.^\gok + W)); for i = j, 



where 



Q = ^°(5oo + W)- (3.5c) 
Substituting Eqs. (3.5) in Eq. (3.4) gives 

CjjF = U^C^F = 0, (3.6a) 

or 

CyF{x'',v') = 0, (3.6b) 

where 

+ E"'(»o* +»'))^}, (3.6c) 

We caution that the post-Newtonian hydrodynamics is obtained from inte- 
grations of Eq. (3.6a) over the v-space rather than Eq. (3.6b) (see appendix 
C). Next we expand jCy up to order (v/c)^- For this purpose, we need 
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expansions of Einstein's field equations, the metric tensor, and the affine 
connections up to various orders. Einstein's field equation with harmonic 
coordinate conditions, g^^T^^, = 0, yields (see appendix A): 

SttG 



y2 2 



y2 4 



S'OO 



OjiOO 



c 

d'^'^goo , 2 d'^'^goo (d'^gm d'^goQ 



2 dt^ 
SttG. 



dx^dx^ 



dx^ 



dx^ 



OrpOO 



- ~r~ ^ ' 



y2 2 



9ij = 



c 

SttG 



,4 *J 



(3.7a) 

(3.7b) 
(3.7c) 
(3.7d) 



The symbols ""g^v and ^T^'^ denote the nth order terms in v/c in the metric 
and in the energy-momentum tensors, respectively. Solutions of Eqs. (3.7) 
are 



'goo = -24>lc\ 

3„ t /„3 



gio = 6/c , 



S'oo 



-2( 



+ ^)/c^ 



(3.8a) 
(3.8b) 
(3.8c) 
(3.8d) 



where 



0(x,i) = 

r(x,t) = 
V'(x,i) = 



G r0T00(x',t)^3,. 



J |x — X 



X — X' 



X — X' 



1 920(x',t) 



47r 



+ G'^T""(x',t) 



+G V*(x',t) 



(3.9a) 
(3.9b) 



(3.9d) 
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where a bold character denotes a three-vector. Substituting Eqs. (3.8) and 
(3.9) in (3.6c) gives 



pn 



d_ 
di 



. d 



d(f) d 



dx'- dx^ dv^ 



1 r/. , 2^94> d(f) ,■ jd(f) dip 

\(4(h + V )— —v v^ - V — H — 



+ 



c^^' ' dx^ dx^ dt dx^ 



dt ^ dv' 



(3.10) 



' dx^ dx^ 

where and are the classical Liouville operator and its post-Newtonian 
correction, respectively. Equation (3.6b) for the distribution function F{x^, 
becomes 

{C'^ + Cf"')F{t,x\v') =0. (3.11) 

The classical Liouville's equation and its symmetries have been studied 
extensively by Sobouti |2T|, |2^, ^ ; Sobouti and Samimi |^ ; Samimi 

and Sobouti Sobouti and Dehghani ||23|; Dehghani and Sobouti p^, ^ . 



The three scalar and vector potentials 0, ip and ^ can now be given in terms 
of the distribution function. The energy-momentum tensor in terms of 
Fix", U') is 



T'^'ix^ 



F{x\U')y/^£U, 



where g = det{g^i,). For various orders of T^'^ one finds 



°T°°(a;^) =c^y F{x\v')d'v, 
2t0°(x^) = J{v^ + 20(x^))F(x\ v')dh, 
^T^{x^) = J vVF{x\v')dh, 
iT°*(x^) = cj v'F{x\v')dh. 
Substituting Eqs. (3.13) in (3.9) gives 



(3.12) 

(3.13a) 
(3.13b) 
(3.13c) 
(3.13d) 

(3.14a) 



27 



C(x, t) = -AG I J_ ^ dV (3.14b) 



.2G / -"^<-''--'),ir- 



+2G 



|x — x'l 

Frx',t,v')Frx",t,v"^ 



/^Af^^qfA^^dr'dr". (3.i4c) 

y X — x' x' — x' 



where dV — d^xd^v. Equations (3.11) and (3.14) complete the pn order 
of Liouville's equation for self gravitating systems embeded in a flat space- 
time. 



3.2 Post-Newtonian Liouville's equation: Static 
solutions 

In the last section we obtained Liouville's equation in pn approximation. In 
this section we seek static solutions of pn/, -F(x, v). In this time-independent 
regime macroscopic velocities along with the vector potential ^ vanish. 
Equations (3.10) and (3.11) reduce to 

One easily verifies that the following, a generalization of the classical energy 
integral, is a solution of Eq. (3.15) 

= l^;2 + 0+(202 + ^)/c2. (3.16) 

The first two terms is exactly classical energy integral and the other terms 
come out from pn correction. Furthermore, if 0(x) and '0(x) are spherically 
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symmetric, which actually is the case for an isolated nomotating system in 
an asymptotically flat space-time, the following generalization of angular 
momenta are also integrals of Eq. (3.15) 

h = SijkX^v'' exp{-(p/c^) ^ EijkX^v^il - ^/c^), (3.17) 

where eijk is the Levi-Cevita symbol. Static distribution functions maybe 
constructed as functions of E and even functions of h. The reason for 
restriction to even functions of /f" is to ensure the vanishing of the 
condition for validity of Eq. (3.15). 

3.3 Post-Newtonian polytropes 

In addition to hydrodynamics equations, one need an equation of state to 
determine comletely a theoretical model for a star. The equation of state 
describe a relation between the mass density and the pressure of the system. 
In order of choosing equation of state, the theoretical models are different 
for a star. Polytropic model is a simple theoretical model to describe the 
equilibrium of star. It relates the pressure to the mass density to power of 
r, the adiabatic index. Classical polytropic model are studied by Eddington 

As in classical polytropes we consider the distribution function for a 
polytrope of index n as 

Fn{E) = -^(-E)"-3/2; for E<0, 

= for E > 0, (3.18) 

where a„ is a constant. By Eqs. (3.13) the corresponding orders of the 
energy-momentum tensor are 

V = anl3nc\-Ur, (3.19a) 
29 



where 



2T°° = 2anM{-Ur + '2anln{-UY^\ (3.19b) 

= = 2a„7„(-t/)'^+\ (3.19c) 

^T^' = 0, (3.19d) 

(5n =/'(!- Cr-'"\"^d(^ = r(3/2)r(n - l/2)/r(n + 1) , (3.20) 
io 

7n =/'(!- C)"~'/'C'/'rfC = r(5/2)r(n - l/2)/r(n + 2) , (3.21) 





and U — (j)-\-2(jP' /(? is the gravitational potential in pn order. It will 

be chosen zero at the surface of the stellar configuration. With this choice, 
the escape velocity — \f—2U will mean escape to the boundary of the 

system rather than to infinity. Einstein's equations, Eqs. (3.7), (3.8) and 
(3.9), lead to 

= ^ Or°° = AT:Ganl3n{-UT, (3.22) 
= 47rG'('T°° T") = ^nGanMi-UY 

+167rGa„7n(-[/)"+^ (3.23) 

Expanding {—UY as 

(-[/)" = {-<t>ni + 71(20 + |)/c^], (3.24) 

and substituting it in Eqs. (3.22) and (3.23) gives 

VV = AT:GaM-4>T - 2n(-0)"+Vc2 - n(-0)"-V/c'] , (3.25) 

V^V' = 47rG'«„/3„(4^ - 2)(-0)"+^ (3.26) 

For further reduction we introduce the dimesionless quantities 

X = a C, (3.27a) 

-0(x) = A^^(C), (3.27b) 

-i^{x) = A^e(C), (3.27c) 

-Cix) = X'^WiO, (3.27d) 
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where, in terms of pc, the central density, A — {Pc/ Cinf^nY^^ and a — 
AnGpc/X. Equations (3.25) and (3.26) reduce to 

V^^ + r = ?n(2^"+^ - ^-^e) , (3.28a) 
V^e + (4^ - 2)^"+^ = 0, (3.28b) 

Pn 

where = ^^(C^^)- dimensionless pn expansion parameter q 

emerges as 

^- -i?2Cii^'(Ci)r ^^-'^^ 

where it!s is the Schwarzschild radius, it! = a^i is the radius of system, and 
Ci is the first zero of ^(C), ^(Ci) — 0- The order of magnitude of q varies 
from 10~^ for white dwarfs to 10~^ for neutron stars. For future reference, 
let us also note that 

- U ^ X[e + q{e - 29'^)]. (3.30) 

We use a forth-order Runge-Kutta method to find numerical solutions of 
the two coupled nonlinear differential Eqs. (3.28). At the center we adopt 

9(0)^1; ^'(0) = ^|o=0. (3.31) 

In tables 3.1 and 3.2, we summarize the numerical results for the Newtonian 
and post-Newtonian polytropes for different polytropic indices and q values. 
The pn corrections tend to reduce the radius of the polytrope. The larger 
the polytropic index and/or q the larger this reduction. 
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Table 3.1: A comparison of the Newtonian and post-Newtonian polytropes 
at certain selected radii for n=l, 2, 3, 4 and 5, and different values of q. 



n 


Polytropic 


Newtonian 


pn polytrope, 9 + g(G — 29'^) 




radius, C, 


polytrope, 9 


q = IQ-^ 


q = 10-^ 


q = 10-^ 




0.0000000 


1.00000 


1.00000 


1.00000 


1.00000 




1.0000000 


0.84147 


0.84147 


0.84156 


0.85043 




2.0000000 


0.45465 


0.45465 


0.45470 


0.46069 


1 


3.0383400 


0.03393 


0.03392 


0.03358 


0.00000 




3.1403800 


0.00039 


0.00038 


0.00000 






3.1415800 


0.00001 


0.00000 








3.1415930 


0.00000 










0.0000000 


1.00000 


1.00000 


1.00000 


1.00000 




2.0000000 


0.52984 


0.52984 


0.53005 


0.55904 




4.0000000 


0.04885 


0.04884 


0.04858 


0.02500 


2 


4.1451500 


0.02776 


0.02775 


0.02746 


0.00000 




4.3501500 


0.00035 


0.00033 


0.00000 






4.3528000 


0.00001 


0.00000 








4.3529000 


0.00000 










0.0000000 


1.00000 


1.00000 


1.00000 


1.00000 




3.0000000 


0.38286 


0.38286 


0.38315 


0.41848 




6.0000000 


0.04374 


0.04373 


0.04338 


0.01817 


3 


6.2838000 


0.02854 


0.02853 


0.02816 


0.00000 




6.8862000 


0.00044 


0.00043 


0.00000 






6.8964000 


0.00001 


0.00000 








6.8967000 


0.00000 
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Table 3.2: Same as Table 3.1 



n 


rolytropic 


AT J_ 

Mewtonian 


pn polytrope, 9 + g(0 — 29^) 




radius, ^ 


polytrope, u 


q = 10^5 


g = 10 


q = 10 ^ 




0.0000000 


1.00000 


1.00000 


1.00000 


1.00000 




3.0000000 


0.44005 


0.44005 


0.44022 


0.46949 




6.0000000 


0.17838 


0.17838 


0.17818 


0.17746 




9.0000000 


0.07955 


0.07954 


0.07919 


0.06496 


4 


12.5013000 


0.02350 


0.02349 


0.02304 


0.00000 




14.0000000 


0.00802 


0.00801 


0.00753 






14.8625000 


0.00051 


0.00050 


0.00000 






14.9705000 


0.00001 


0.00000 








14.9713400 


0.00000 










0.0000000 


1.00000 


1.00000 


1.00000 


1.00000 




5.0000000 


0.28480 


0.28480 


0.28482 


0.29394 




10.0000000 


0.11894 


0.11894 


0.11862 


0.10940 


4.5 


12.2000000 


0.08779 


0.08779 


0.08743 


0.00000 




15.0000000 


0.06125 


0.06125 


0.06085 






20.0000000 


0.03231 


0.03230 


0.03185 






25.0000000 


0.01498 


0.01492 


0.01444 






30.0000000 


0.00334 


0.00333 


0.00284 






31.2256000 


0.00107 


0.00106 


0.00000 






31.7847000 


0.00001 


0.00000 








31.7878400 


0.00000 
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Chapter 4 

The post Newtonian modes 



In the last chapter we obtained Liouville's equation in pn approximation. 
Furthermore, we found the integrals of pnl, generalization of the classical 
energy and angular momentum, and constructed an equilibrium distribution 
function for the system. In this chapter, we study the non-equilibrium state 
of a stellar system in pn approximation. We assume a small perturbation 
in the system, i.e. in the distribution function, and obtain the linearized 
Liouville's equation. Finally, using the linearized equation, we study normal 
modes of the system in pn approximation. In this chapter all quantities are 
dimensionless. 



In section ^]T] we give the pn order of the linearized Liouville equation 
that governs the evolution of small perturbations from an equilibrium state. 



In sections [4.2| and we extract the equation for a sequence of new modes 
that are generated solely by pn force but are absent in classical regime. In 
section O we explore the 0(3) symmetry of the modes and classify them 



on basis of this symmetry. In section ^]5| we study hydrodynamics of these 



modes. In section we seek a variational approach to the calculation of 
pn modes and give numerical values for polytropes. 
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4.1 Post Newtonian Liouville'e equation, Lin- 
earized 



In chapter 3 we obtained Liouville's equation in the post-Newtonian ap- 
proximation (pnl) for the one particle distribution of a gas of coUisionless 
particles as 

(-z^ + £)F(x, u, t) = {-i^ + C^' + gr")F(x, u, t) = 0, (4.1) 

where (x, u) are phase space coordinates, g is a small post-Newtonian ex- 
pansion parameter, the ratio of Schwarzchild radius to a typical spatial 
dimension of the system, Eq. (3.29). The classical and post-Newtonian 
operators, C'^'' and respectively, are 

+ (4-2a) 

(4.2b) 

The imaginary factor i is included for later convenience. The potentials 
9{'x,t), 6(x, t) and /^(x, t), solutions of Einstein's equations in pn approxi- 
mation, are 

F(x',t,u')^^, , , ru'F(x',t,u') 



e(x,t) = / ; ' ' / dV, 77 x,t =4 / — / ' U T', 4.3a, b 
J |x — x'l J |x — x'l 

^ ' ^ An J |x-x'||x'-x"| J |x-x'| 

f F(x.'.t.u')F(x.".t.u") , „ , , 

-2 / \ ' ^ „ ^ d^dr, 4.3c 

J |x — x'||x — x"| 

where ciF = d^xd^u. See chapter 3 for details. In an equilibrium state, 
F(x, u) is time-independent. If, further, it is isotropic in u, macroscopic 
velocities along with the vector potential 77 vanish. It is also shown in chap- 
ter 3 that the following generalizations of the classical energy and classical 
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angular momentum are integrals of pnl: 



k = eijkX^u''exp{qe) « lf{l + q9), 



(4.4a) 
(4.4b) 



for spherically symmetric 9{r) and ©(r). Equilibrium distribution functions 
in pn approximation can be constructed as appropriate functions of these 
integrals. In chapter 3 the pn models of polytrope were studied in this 
spirit. 

Here we are interested in the time evolution of small deviations from 
a static solution. Let F F{e) +5-F(x, u, t), | 5F F, V x, u, t. 
Accordingly, the potentials split into large and small components, 9{r) + 
59{x,t), 6(r) + 5Q{x,t) and 5r]{x,t) where r = |x|. Both the large and 
small components, can be read out from Eqs. (4.3). Substituting this 
splitting in Eq. (4.1) and keeping terms linear in 5F gives 



where C is now calculated from Eqs. (4.2) with 9{r), 0(r) and 77 = 0. 



i—SF = CSF + SCF{e), 



(4.5) 



Thus 



£ = 



+ qjc:' 



(4.6a) 





(4.6b) 




d 



(4.6c) 



For SC Eqs. (4.2), similarly, give 




(4.7a) 
(4.7b) 



2 d69 j dSrji 



(4.7c) 



dt dt 
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Equations (4.5)-(4.7) are the generalizations of the hnearized classical Liouville- 
Poisson equations to pn order. The classical case was studied briefly by 
Antonov p^. He separated 6F into even and odd components in u and 



extracted an eigenvalue equation for 5Fodd- Sobouti ||2T], |2^, ^ ^ elab- 
orated on this eigenvalue problem, studied some of its symmetries and ap- 
proaches to its solution. Sobouti and Samimi and Samimi and Sobouti 
showed that Antonov's equation has an 0(3) symmetry and its oscilla- 



tion modes can be classifled by a pair of eigennumbers (j, m) of a pair phase 
space angular momentum operators (J^, J^). In analyzing Eqs. (4.5)-(4.7) 
we have heavily relied on these studies. 



4.2 The Hilbert space 

Let Ti be the space of complex square integrable functions of phase coordi- 
nates (x, u) that vanish at the phase space boundary of the system: 

n : /(x, u); j f*f^/^dV = flnite, /(boundary) = 0, (4.8) 

where y/—g = 1 + 2q6 in pn order. Integrations in Ti are over the volume 
of the phase space available to the system. In particular the boundedness 
of the system sets the upper limit of u at the escape velocity V^, where 
= ^(x) is the gravitational potential at x. Thus, /(x, ^J29{x)) = 0. 
Theorem : C = C^^ + g£^" of Eqs. (4.6) is Hermitian in 7i, 



I g\lif) (1 + 2qd)dY = J {Cg)*/ (1 + 2qe)dT; gJeH (4.9) 

Proof: Substituting Eqs. (4.6) in (4.9), carrying out some integrations by 
parts over the x and u coordinates and letting the integrated parts vanish 



37 



on the pahse space boundary: 

I g*{Cf) (1 + 2qe)dT = j g*{C'' + qC^f (1 + 2g^)rfr, 

= j g*c^^fdr + q(2 j e g*C^^fdr + J g*CP" fdr^ . 

(4.9a) 

At the classical order, the classical Liouville operator, C'^'-, is Hermitian in 



n, [21 



J g*{C^'f) dV = j {C^'gYf dV- gJeU 

Therefore the first two terms in RHS of Eq. (4.9a) will 
be 

j{C'grfdT + 2qj[C''' mrfdV, 

= J {C^gYfdV + 2qJ {C^gYfOdV + 2tqJ^-u e'g*fdT. (4.9b) 

The third term, the post-Newtonian Liouville operator, at the pn order 
is not Hermitian, then 

j g*CP" fdr = J {D^'gYfdT - 2tq J ^ ■ u e'g*fdT (4.9c) 

The proof will be completed by adding Eqs. (4.9b) and (4.9c), QED. 

The term 5C is not, in general, Hermitian. Nonetheless, one may proceed 
as Antonov did with the classical case and obtain a second order differential 
operator (almost square of C + 5C) in some subspace of Ti. We are, however, 
pursuing a much simpler problem here in which 5C term vanishes identically 
leaving Eq. (4.5) as an eigenvalue problem governed with the Hermitian 
operator C alone. 



4.3 The post-Newtonian modes 

The effect of pn corrections on the classical solutions of Eq. (4.5) can be 
analyzed by the usual perturbation techniques. Whatever the procedure, 
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the first order corrections on the known modes will be small and will not 
change their nature. We will not pursue such issues here. The main interest 
of this work is to study a new class of solutions of Eq. (4.5) that originate 
solely from the pn terms and have no precedence in classical theories. It 
is not difficult to anticipate the existence of such modes. Perturbations on 
an equilibrium state, that are functions of classical integrals (energy and 
angular momentum, say) do not disturb the equilibrium of the system at 
classical level. That is they do not induce restoring forces in the system. 
They, however, do so in the pn regime, and make the system oscillate about 
the pn equilibrium state. Such perturbations may be considered as a class 
of infinitely degenerate zero frequency modes of the classical system. The 
pn forces unfold this degeneracy and turn them into a sequence of non zero 
frequency modes distinct and uncoupled from the other classical modes. We 
have termed them as pn modes. 

A hydrodynamic interpretation of pn modes is the following. In spher- 
ically symmetric fiuids, toroidal motions are neutral. Sliding one spherical 
shell of fiuid over the other is not opposed by a restoring force. The pn 
forces or for that matter a small magnetic field or a slow rotation (mainly 
through Coriolis forces) gives rigidity to the system. The fiuid resists against 
such displacements and a sequence of well defined toroidal modes of oscilla- 
tion develop. See Sobouti Hasan and Sobouti [^, Nasiri and Sobouti 



45[] , and Nasiri |^ for examples and typical calculations in the case weak 



magnetic fields and slow rotations. 

In the Fourier time transform of Eq. (4.5), 

£6F + 6CF{e) = u5F, (4.10a) 

we split 5F into even and odd terms in u. Thus, 

5F(x,u) =G'_(x,u) + G'+(x,u), G'±(x, u) = ±G±(x, ±u). (4.10b) 
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Considering the fact that both L and are odd in u, Eq. (4.10a) sphts 
accordingly: 

LG- + qujFeU^Se = u;G+, (4.11a) 
jCG+ - iF^u'-^ [56 + q{5Q - mO)] - qujF^u'5'qi = tuG_, (4.11b) 

where 

G^(xf,n')^^, _ ^ /-u'G^x'.u') 



y |x — x'l y |x — x'l 

^ ' ^ 47ry |x-x'||x'-x"| i 

-2 / G.(x-,u')F(e") + F(e')G,(x",u")^p^p, 
J |x — x'l |x' — x"| 

Operating on Eq. (4.11a) by C and substituting for jCG+ from Eq. (4.11b) 

gives a second order differential equation for G_: 

C^G^ = u^G^+iuF^u'^ [59 + q(5e - Ae5e)\+quj^F^u'5r]i-qujFMv?5e). 

(4.13a) 

We now seek a solution of Eq. (4.13a) in the form of classical energy and 
angular momentum integrals, G'_(x, u) = G_(e"^',Zf). In the next section, 
after we discuss the 0(3) of Eq. (4.13a), we show that such solutions can 
be chosen from among the eigenfunctions of a pair of phase space angular 
momentum operators, (J^, J^). We also show that for such solutions 56 and 
50 vanish identically reducing Eq. (4.13a) to 

jC^G^ = cu^ (G. + qFeu'5r]i) . (4.13b) 

Multiplying Eq. (4.13b) by G*_, integrating over the phase space volume of 
the system, and considering the facts that C — C'^ + qC^ is Hermitian and 
C'^G_{e''\lf) = 0, gives 

j {CG^yCG^l + 2q6)dV ^q^ J {a"'GJ)*a"'G^{l + 2q6)dV 
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^2 



J G*_G4l + 2qe)dT + qj G*_F,u'67]i{l + 2qe)dV . 

(4.14a) 



j {Cf"'G^yCf''G^dT = ^ J G*_G^dT. (4.14b) 



Equation (4.14a) shows that uj is of the same order of smallness as q. Thus, 
ehminating the terms of order q^, uj'^q and higher reduces Eq. (4.14a) to 

Equation (4.14b) provides a variational expression for uj"^ and will be used 
as such to calculate the allowable uo"^. The frequencies, uo, are real mean- 
ing that the corresponding deviations from the equilibrium state are stable 
oscillation modes. Furthermore, these perturbations will be different from 
the conventional classical modes, for they are excited by pn terms in the 
equations of motion that are absent at classical level. 

4.4 0(3) symmetry of £ = C^^ + qC^^ 

For spherically symmetric potentials, 9{r) and 0(r), both C^^ and de- 
pend on the angle between x and u and their magnitudes. Simultaneous 
rotations of the x and u coordinates about the same axis by the same angle 
leaves these operators form invariant. The generator of such simultaneous 
infinitesimal rotations on the function space Ti is 

J, = 4 = -ie,,k [x^^ + n^'^j , (4.15) 

which has the angular momentum algebra 

[Ji, Jj] = iSijkJk- (4.16) 



Commutation of Jj with was first established by Sobouti pT|. Here we 



confine the discussion to the symmetry of O'^. Straightforward calculations 
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reveal that 

[/:^",J.] = 0, (4.17) 



since 



r ou'^ 



Thus, it is possible to choose the eigensolutions, G_ of Eq. (4.14b) si- 
multaneously with those of and Jz- The eigensolutions of the latter 
pair of operators are worked out in the appendix D. They are of the form 



/(e^', lf)Ajm', i^rn integers, where / is an arbitrary function of the classical 
integrals and Kjm is a complex poljTiomial of order j of the components of 
the classical angular momentum. If . The x and u parity of Kjm is that of 
j. See appendix D for proofs this statement. 

We are now in a position to point out an interesting feature of the 
eigenmodes. Both uj"^ and in Eq. (4.13b) and the integrals in Eq. (4.14b) 
are real. Thus, G_ can be chosen real or purely imaginary. By Eq. (4.11a), 
then will be purely imaginary or real. That is, an eigensolution 5F = 
G-+G+ belonging to a nonzero u; is a complex function of phase coordinates 
in which both the x and u parities of the real and imaginary parts are 
opposite to each other. This feature is shared by the classical modes of the 
classical Liouville's and Antonov's equation. 

In section we will take a variational approach to solutions of Eq. 



(4.14b). As variational trial functions we will consider the following 

N 

G- = fjm = f{e)Re Ajm = { c„(-e)"]i?eAj>„, j = odd, c„ = consts. 

n=j+l 

(4.18) 
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Combining this with its corresponding even counterpart from Eq. (4.10a) 
we obtain 

SF,^{^, u, t) = {l + ^Cnfjme-'^'- (4.19) 

At this stage let us note an important property of Liouville's equation. If a 
pair [u, 5F) is an eigensolution of Liouville's equation, {—u, 5F*) is another 
eigensolution. This can be verified by taking the complex conjugate of 
Eq. (4.10a). These solutions, being complex quantities, cannot serve as 
physically meaningful distribution functions. Their real or imaginary parts, 
however, can. With no loss of generality we will adopt the real part. Thus, 

Re 5Fjm{^, u, t) = f{e)Re Aj^ cosut + i-a"'{f{e)Re Ajm) sinuot. (4.20) 

The eigenmodes of Eq. (4.10a) are m-independent. By m-independence 
we mean a) the eigenvalues cu do not depend on m and are 2j + 1 fold 
degenerate, and b) the expansion coefficients, c„, of Eq. (4.12) do not 
depend on m. Proof: From the appendix D, Eq. (D.4), J± = ± iJy 
are ladder operators for Ajm- Operating on fjm of Eq. (4.18) by J± will 
give the mode /j,m±i without changing the expansion coefficients. Secondly, 
substituting J±fjm — \J (j ^ m) (j ± m + l)/j,m±i in Eq. (4.14a) instead of 
fjm, and noting that fjmS can be normalized for all m's, cu^ will remain 
unchanged. 

4.5 Hydrodynamics of pn modes 

In this section we calculate the density fluctuations, macroscopic velocities, 
and the perturbations in the space-time metric generated by a pn mode. It 
was pointed out earher that for j an odd integer, fjm{^, u) of Eq. (4.18) is 
odd while jC^fj^ is even in both x and u. The macroscopic velocities are 
obtained by multiplying Eq. (4.20) by u and integrating over the u-space. 

43 



Only the odd component of SFjm contributes to this bulk motion, 

pv = j f{e) Re A j^ud^u cos ujt. (4-21) 

In appendix D, Eqs. (D.ll), we show that pv is a toroidal spherical har- 
monic vector field. In spherical polar coordinates it has the following form 

—1 d dYm 
p{vr, v^, v^) ^r^G{ves){0, Re - — r— Y,>n(t?,<p), Re cosujt, 

sm if 0(p ov 

(4.22a) 

where 

G{ves)= r f{e)u^^^du, (4.22b) 
Jo 

and Ves — V26 is the escape velocity from the potential 9(r). The macro- 
scopic density, generated by the even component of Eq. (4.20), is 

Sp{x,t) = i- / a"'{f{e)Re Ajjd^usmcut 

UJ J 

= • J f{e)Re Aj^nd^ smut = 0. (4.23) 

The second integral is obtained by an integration by parts. The vanishing 
of it comes about because of the fact that the radial vector x is orthogonal 
to the toroidal vector pv. One also notes that V • (pv) = 0. It can further 
be verified that, the continuity equation is satisfied at both classical and pn 
level. 

To complete the reduction of Eqs. (4.13) we should also show that 69 
and 50 vanish. The former is zero because 5p — 0. For the latter, from Eq. 
(4.3c) and Eq. (4.20) for SF, one has 

An J |x — x'l J |x — x'l |x' — x"| 

+2 J J u"5F(x', n')d\' = 0. (4.24) 

The vanishing of the first two terms is obvious. The third term vanishes 
because the integral over u' has the same form as in 5p except for the 
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additional scalar factor u'^. Like 5p it can be reduced to the inner product 
of the radial vector x and a toroidal vector. QED. 

The toroidal motion described here slides one spherical shell of the fluid 
over the other without perturbing the density, the Newtonian gravitational 
field and, therefore, the hydrostatic equilibrium of the classical fluid. In 
doing so, it does not affect and is not affected by the conventional classical 
modes of the fluid at this first pn order. 

Nonetheless, the pn modes are associated with space time perturbations. 
Prom Eq. (3.8c) and Eq. (4.3b), goi component of the metric tensor is 

goi^Vi^^f P^rd'x'. (4.25) 
J |x — x'l 

In spherical polar coordinates, one obtains 

Tjr = 0, (4.26a) 
1 d 

rj^ = -UjRe - — ^Tr-^jm{i^^ cosut, (4.26b) 
smiydip 

Tjip = ajRe (i), ip) cos out, (4.26c) 

where 
167r 



ir/Ryyj{R) + (2j + 1^ jj^ r'-^-'yjir')dr' for r < i? 
2j + 1 (R/ry+%(R) iorr> R 

y.{r) = r--'-' r r'^'^^G{9{r'))dr', (4.26e) 

•J 

Gieir)) = / f(eW+^du 
Jo 

= 2^'/2+ir(j/2 + 2)r(n + l)^(r)"+^/2+2/r(n + j/2 + 3), (4.26f) 

where R is the radius of the system and r(n) is the gamma function. The 
remaining components of the metric tensor remain unperturbed. 
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4.6 Variational solutions of pn modes 



We substitute the trial function of Eq. (4.18) in Eq. (4.14b) and turn it 
into a matrix equation. Thus 

^2 

C^WC = —C^SC, (4.27) 

where C = [c„] is the column matrix of the variational coefficients of Eq. 
(4.18), and the elements of S and W matrices are 

Spg = j {-ef+'i\ReKjm?dV, (4.28a) 
Wp, = J (£^"(-e)^i?e A^.^nC^'^i-eyRe Ajn,)dr. (4.28b) 

Minimizing cu^ with respect to variations of C gives the following matrix 
equation 

WC = ^SC. (4.29) 
Eigen w's are the roots of the characteristic equation 

^2 

\W-—S\=Q. (4.30) 

For each u, Eq. (4.29) can then be solved for the eigenvector C. This 
completes the Rayleigh-Ritz variational formalism of solving Eq. (4.14a). 
In what follows we present some numerical values for polytropes. 

4.6.1 pn Modes of polytropes belonging to (j, m) = 

(l,m) 

We analyze the case m = 0, only. From the m-independence of eigen- 
modes (see theorem of section [4.4|) the eigenvalue and the expansion coef- 
ficients, Cn, for m = ±1 will be the same. From Eqs. (D.9), Aiq = 1^ = 
rasin'i?sinQ;sin(/? — ip), where ("i?, v^) and (a,/3) are the polar angles of x, 
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of u, respectively. Substituting this in Eqs. (4.28) and integrating over 
directions of x and u vectors and over < u < V26 gives 

Spg= r Vdx, (4.31a) 



Wpg = TTGpc{{l6apg - bpg) ^ 9'^ eP^'^^^-^ dx 

J 

+(1 - %apg) f e'e'e^+'^+^-^xUx 

J 

+apq C e'^e^^'i+^-^x^dx}, (4.31b) 
Jo 

^ pq{p + q + 2.5) 

(4.31c) 

_ 4(p + g)^ + 9(p + g)-13 ^ ^ ^ 
(p + ?-l)(p + g + 3.5) ^ 

Polytropic potentials 6 and © were obtained from integrations of Lane Em- 
den equation and Eqs. (3.28), respectively. Eventually, the matrix elements 
of Eqs. (4.31), the characteristic Eq. (4.30) and the eigenvalue Eq. (4.29) 
were numerically solved in succession. Tables 4.1-4.4 show some sample 
calculations for polytropes 2, 3, 4, and 4.9. Eigenvalues are displayed in 
lines marked by an asterisks. The column following an eigenvalue is the 
corresponding eigenvector, i. e. the values of Ci, C2, • • •, of Eq. (4.18). To 
demonstrate the accuracy of the procedure, calculations with six and seven 
variational parameter are given for comparison. The first three eigenval- 
ues can be trusted up to two to four figures. Convergence improves as the 
polytropic index, i.e. the central condensation, increases. Eigenvalues are 
in units of irGpcq^ and increase as the mode order increases. 
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Table 4.1: pn modes of polytrope n=2, belonging to (j, m) = (1, 0). Eigen- 
values are in units AirGpcq^, c„'s are the linear variational parameters of 
Eq. (4.22). A number a x 10^^ is written a a ± 6. To appraise the accu- 
racy of the computations two sets of data with six and seven variational 
parameters are given. The first three eigenvalues are reliable up to three 
figures. Characteristically, the accuracy deteriorates as one goes to higher 
order modes. 





.1825+01 


.4973+01 


.6448+01 


.1216+02 


.3425+02 


.1686+03 




Cl 


.3113+02 


-.8912+02 


.1663+03 


.1344+03 


.7545+01 


-.1399+04 




C2 


.3908+02 


.1045+04 


-.3234+04 


-.9746+03 


-.2392+04 


.8484+04 




C3 


-.1420+03 


-.6649+04 


.1801+05 


.4514+04 


.7952+04 


-.9647+04 




C4 


.5803+03 


.1804+05 


-.4351+05 


-.7014+04 


-.2607+03 


-.2251+05 




C5 


-.9110+03 


-.2210+05 


.4724+05 


.8324+03 


-.1811+05 


.5188+05 




C6 


.5252+03 


.1020+05 


-.1874+05 


.2882+04 


.1317+05 


-.2717+05 






.1823+01 


.4865+01 


.5895+01 


.9113+01 


.1465+02 


.4228+02 


.3226+03 


Cl 


.3028+02 


-.7086+02 


.1529+03 


-.3129+02 


.1561+03 


-.4624+02 


.2042+04 


C2 


.4812+02 


.6908+03 


-.2810+04 


.1313+04 


-.1513+04 


-.2762+04 


-.1461+05 


C3 


-.1305+03 


-.3993+04 


.1702+05 


-.5686+04 


.6685+04 


.1077+05 


.2271+05 


C4 


.2576+03 


.8181+04 


-.4788+05 


.3425+04 


-.3673+04 


.1875+04 


.4154+05 


C5 


.1303+03 


-.3086+04 


.6823+05 


.2433+05 


-.2910+05 


-.4718+05 


-.1496+06 


C6 


-.7534+03 


-.7924+04 


-.4771+05 


-.4855+05 


.5132+05 


.5873+05 


.1425+06 


C7 


.5475+03 


.6707+04 


.1302+05 


.2568+05 


-.2386+05 


-.2120+05 


-.4423+05 



prii pn2 pn^ pn^ pn^ pn^ prij 
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Table 4.2: Same as Table 4.1. n = 3 and (j, m) = (1, 0). 





.1534+01 


.4836+01 


.9473+01 


.1938+02 


.4083+02 


.1128+03 




Cl 


.9752+02 


-.6975+02 


.2464+03 


-.2246+03 


-.9102+03 


.3169+04 




C2 


.3284+02 


-.8725+03 


-.1121+04 


-.2590+04 


.1713+05 


-.2631+05 




C3 


.2096+03 


.3859+04 


.5591+04 


.1444+05 


-.1023+06 


.6390+05 




C4 


-.5354+03 


-.5728+04 


-.1216+05 


-.9903+04 


.2599+06 


-.3406+05 




C5 


.3941+03 


.2528+04 


.5215+04 


-.2221+05 


-.2933+06 


-.4814+05 




C6 


.1803+01 


.1125+04 


.3307+04 


.2153+05 


.1208+06 


.4268+05 




U? 


.1533+01 


.4688+01 


.7993+01 


.9068+01 


.1124+02 


.1909+02 


.1093+03 


Cl 


.9318+02 


-.1440+03 


-.1202+03 


-.1069+04 


-.5706+03 


-.5482+02 


.3703+04 


C2 


.1121+03 


.6997+03 


.5482+04 


.1856+05 


.7685+04 


-.5626+04 


-.3381+05 


C3 


-.2118+03 


-.4506+04 


-.2955+05 


-.1063+06 


-.4112+05 


.3078+05 


.1007+06 


C4 


.2709+03 


.9777+04 


.5298+05 


.2726+06 


.7791+05 


-.4371+05 


-.1109+06 


C5 


.1206+03 


-.9309+03 


-.6283+04 


-.3375+06 


-.1278+05 


-.7049+03 


.1239+05 


C6 


-.7005+03 


-.1574+05 


-.7154+05 


.1894+06 


-.9027+05 


.3228+05 


.4581+05 


C7 


.5309+03 


.1200+05 


.5087+05 


-.3511+05 


.5945+05 


-.1218+05 


-.1722+05 



prii pris ^"-4 V^h P^e pn-j 
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Table 4.3: Same as Table 4.1. n = 4 and (j, m) = (1, 0). 





.7569+00 


.2822+01 


.5661+01 


.8814+01 


.1519+02 


.6952+02 




Cl 


.6291+03 


-.1067+04 


.2143+04 


-.1949+04 


-.6870+04 


.1400+05 




C2 


-.9217+02 


.1770+04 


-.1693+05 


.1131+05 


.8373+05 


-.2337+06 




C3 


.4162+03 


.2808+04 


.5682+05 


-.3654+04 


-.3195+06 


.1293+07 




C4 


-.3883+04 


.5860+04 


-.1184+06 


-.2807+05 


.4791+06 


-.3112+07 




C5 


.6427+04 


-.2303+05 


.1257+06 


-.4668+04 


-.2545+06 


.3371+07 




C6 


-.3089+04 


.1612+05 


-.4514+05 


.3416+05 


.1251+05 


-.1344+07 






.7569+00 


.2813+01 


.5021+01 


.8747+01 


.1272+02 


.3322+02 


.7683+02 


Cl 


.5590+03 


-.8716+03 


.2653+03 


-.2421+04 


.1881+04 


.1412+05 


.3376+05 


C2 


.1189+04 


-.2018+04 


.1406+05 


.1926+05 


-.7436+04 


-.2356+06 


-.5191+06 


C3 


-.6377+04 


.2349+05 


-.1057+06 


-.4732+05 


-.5363+05 


.1298+07 


.2528+07 


C4 


.9376+04 


-.3509+05 


.2059+06 


.6165+05 


.2228+06 


-.3112+07 


-.4750+07 


C5 


.5449+03 


-.4645+04 


-.2977+05 


-.4272+05 


-.7106+05 


.3356+07 


.2298+07 


C6 


-.1192+05 


.4364+05 


-.2533+06 


-.3854+05 


-.4046+06 


-.1333+07 


.2455+07 


C7 


.7228+04 


-.2275+05 


.1775+06 


.5845+05 


.3227+06 


-.1382+03 


-.2085+07 



prii pris ^"-4 V^h P^e pn-j 
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Table 4.4: Same as Table 4.1. n = 4.9 and {j, m) = (1, 0). 





.4481+00 


.1827+01 


.4078+01 


.6515+01 


.1170+02 


.1391+03 




Cl 


-.2888+02 


.1663+03 


-.2794+03 


.1593+03 


.1405+03 


.1081+05 




C2 


-.2440+03 


-.7593+04 


.2050+05 


-.2099+05 


.2665+05 


-.2129+06 




C3 


.4933+05 


-.2772+04 


-.1400+06 


.1883+06 


-.3467+06 


.1344+07 




C4 


-.1722+06 


.1443+06 


.2902+06 


-.5138+06 


.1372+07 


-.3583+07 




C5 


.2124+06 


-.2675+06 


-.2194+06 


.4871+06 


-.2092+07 


.4207+07 




C6 


-.8916+05 


.1394+06 


.5712+05 


-.1179+06 


.1073+07 


-.1790+07 






.4380+00 


.1805+01 


.4006+01 


.6190+01 


.7980+01 


.1439+02 


.8964+02 


Cl 


-.1701+02 


.1379+03 


-.3341+03 


.3427+03 


-.3020+03 


.7695+03 


.8642+04 


C2 


-.6649+03 


-.6322+04 


.2326+05 


-.3097+05 


.2196+05 


-.1111+05 


-.1534+06 


C3 


.5135+05 


-.1143+05 


-.1601+06 


.2940+06 


-.2552+06 


.1349+06 


.8174+06 


C4 


-.1667+06 


.1599+06 


.3264+06 


-.9227+06 


.1022+07 


-.9712+06 


-.1565+07 


C5 


.1694+06 


-.2551+06 


-.1784+06 


.1132+07 


-.1574+07 


.3018+07 


.4968+06 


C6 


-.1770+05 


.8656+05 


-.9582+05 


-.4879+06 


.7432+06 


-.3959+07 


.1421+07 


C7 


-.3646+05 


.3341+05 


.9586+05 


.2938+05 


.8318+05 


.1819+07 


-.1047+07 



pni pn2 pris pn^ pn^ puQ pn-j 



51 



Part II 

The stability curve of r-modes 
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Chapter 5 
Introduction 



Rotating neutron stars and black holes have been the objects of many as- 
trophysical studies in recent years. Their strong gravitational fields make 
them ideal laboratories for testing predictions of the theory of general rel- 
ativity. Both types of compact objects are interesting for different reasons. 
Black holes are objects which have completely collapsed under their own 
gravitational field. Since they are curved empty space, black holes are rel- 
atively simple objects to describe. Most observed phenomena from quasars 
(young, active galaxies) can be explained consistently by assuming a cen- 
tral supermassive black hole exists at the galaxy's core. In contrast neutron 
stars, possibly densest configurations of matter which are stable to grav- 
itational collapse, have more complicated structures. Their study against 
requires a diverse range of physics. The interior structure of a neutron star 
includes such features as a normal fiuid coexisting with superfluidity and 
superconductivity, various nuclear processes, rapidly rotating configuration, 
strong magnetic fields, and many other features. They are one of the most 
fascinating objects of theoretical investigation. Observational features, such 
as strong X-ray emission, periodically pulsating radio waves in a relatively 
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narrow beam, sudden spinning up of rotational frequency (glitch), open a 
window for deep understanding of their internal structure. 

Neutron stars are probably among the most promising sources of de- 
tectable gravitational waves by the future generation of gravitational wave 
detectors. Studies suggest that the emitted gravitational radiation by neu- 
tron stars in the Virgo cluster could be detected by the laser interferometer 
gravitational wave detectors such as GEO600, the advanced Laser Interfer- 
ometer Gravitational wave Observatory (LIGO), and VIRGO ^0|. By 2001, 
both LIGO and GEO600 are expected to be sensitive to bursts of amplitude 
around 10~^^. VIRGO, with an even better sensitivity, could come online 



in 2002 or later 61 



Recently Andersson |£2[ numerically showed that in rotating systems 
perturbations driven by Coriolis forces can be unstable at arbitrarily small 
angular velocities. He considered a relativistic but slowly rotating configura- 
tion and found that toroidal perturbations, in the absence of fluid viscosity, 
are unstable because of the emission of gravitational radiation for all rates 
of rotation. These modes which are the relativistic analouge of the New- 
tonaian r-mo des H, H, are unstable even for very slowly rotating perfect 
fluid stars. Thus, the r-mode instability is different to other mode insta- 
bility like /-mode of the star that set at a certain rate of rotation [Q-||7T 



Friedman and Morsink [^] used the relativistic axisymmetric background 
(with slowly rotating approximation) and showed that analytically the in- 
stability is generic: " every r-mode is in principle unstable in every rotating 
star, in the absence of viscosity " . The mechanism of the instability can be 
understood by the generic argument for the gravitational radiation driven 



instability, so-called CFS-instability, after Chandrasekhar [^, Friedman 
and Schutz ||7^, and Friedman |7^, who studied it first. 

Further work which included the effect of viscosity showed the gravi- 
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tational radiation driven instability of the Coriolis modes is important in 
the class of neutron stars which are born with rapid rotations (such as the 
pulsar found in the supernova remnant N157B). 

The excitement over the r-mode instability has generated a large lit- 



erature in the past two years. Andersson et al. and , Lindblom et 
al. independently computed that slowing down of a rapidly rotating, 
newly born star to typical periods of Crab like pulsar 19 ms) can be 
explained by the r-mode instability. This is due to the emission of current- 
quadrupole gravitational radiations, which reduces the angular momentum 
of the star. Kokkotas and Stergioulas |7S| investigated analytically r-mode 
instability for a uniform density Newtonain star and calculated the cor- 
responding timescales and stability curve associated with r-mode. Lind- 



blom, Owen, and Morsink [77| also evaluated the r-mode growing/damping 
timescale by considering fluid viscosity and calculated critical angular ve- 
locities for a polytropic neutron star model. They showed that the coupling 
of gravitational radiation to the r-modes is sufficiently strong to overcome 
internal fluid dissipation effects and so drive these modes unstable in hot 
young neutron stars. This result which has been verified by Andersson, 
Kokkotas, and Schutz |7B[, seemed somewhat surprising at first because 
the dominant coupling of gravitational radiation to the r-modes is through 
the current multipoles rather than the more familiar and usually dominant 
mass multipoles. But it is now generally accepted that gravitational radia- 
tion does drive unstable any hot young neutron star with angular velocity 
greater than about 5% of the maximum (the angular velocity where mass 
shedding occurs). This instability therefore provides a natural explanation 
for the lack of observed very fast pulsars associated with young supernovae 
remnants. Kojima suggested that, in contrast to Newtonian theory, 
r-mode frequencies in general relativity become continuous. This fact is 
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verified mathematically by Beyer and Kokkotas |8( 

The r-mode instability is also interesting as a possible source of grav- 
itational radiation. In the first few minutes after the formation of a hot 
young rapidly rotating neutron star in a supernova, gravitational radiation 
will increase the amplitude of the r-mode (with spherical harmonic index 
m = 2) to levels where non-linear hydrodynamic effects become important 
in determining its subsequent evolution. While the non-linear evolution of 
these modes is not well understood as yet, Owen et al. have developed a 
simple non-linear evolution model to describe it approximately. This model 
predicts that within about one year the neutron star spins down (and cools 
down) to an angular velocity (and temperature) low enough that the in- 
stability is again suppressed by internal fluid dissipation. All of the excess 
angular momentum of the neutron star is radiated away via gravitational 
radiation. Owen et al. |8l[| estimated the detectability of the gravitational 
waves emitted during this spindown, and found that neutron stars spinning 
down in this manner may be detectable by the second-generation ("en- 
hanced") LIGO interferometers out to the Virgo cluster. Bildsten |^ and 
Andersson, Kokkotas, and Stergioulas [^] have raised the possibility that 
the r-mode instability may also operate in older colder neutron stars spun 
up by accretion in low-mass x-ray binaries. The gravitational waves emitted 
by some of these systems (e.g. Sco X-1) may also be detectable by enhanced 
LIGO Thus, the r-modes of rapidly rotating neutron stars have become 
a topic of considerable interest in relativistic astrophysics. 

Furthermore, the r-mode observations open a rich prospect for gravita- 
tional astronomy. Identifying hidden or unnoticed supernovas would be the 
most exciting use of r-mode observation. Another implication of r-mode 
observation is detection of background radiation. If we assume that most 
neutron stars are born with rapid rotation, the background spectrum will 
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reveal not only the star formation rate in the early universe, but also tell us 
about the distribution of initial rotation speeds of neutron stars. Investigat- 
ing r-mode events associated with known supernovas is the other prospect 
of the r-mode observation. This would give us some information about 
cooling rates, viscosity, crust formation, the equation of state of neutron 
matter, and the onset of superfluidity in neutron stars pT| . 

In spite of the recent improvements in our understanding of this insta- 
bility, it seems that the fundamental properties of these modes have not yet 
been sufficiently understood. Previous investigations of the r-modes are re- 
stricted to the case of uniformly and slowly rotating, isentropic, Newtonian 
stars [j76|. A few recent studies were done for relativistic stars with slowly 
rotating and Cowling approximations ||62|. In this sense, it is interesting 
to study the properties of the r-mode instability in the more general cases, 
for example, differentially and rapidly rotating, non-isentropic relativistic 
stars. 

Furthermore, they used a thermodynamic model for the neutron star 
fluid that is not compatible with special relativity, they largely ignored 
superfluid and magnetic field effects. 

In the first work, we address one of r-mode's weaknesses, by utilizing 
a thermodynamic model for the neutron star fiuid that takes the coupling 
between vorticity and shear viscosity into account. Navier-Stokes theory 
has been used to calculate the viscous damping timescales and produce a 
stability curve for r-modes in the T) plane. In Navier-Stokes theory, 
viscosity is independent of vorticity, but kinetic theory predicts a coupling 
of vorticity to the shear viscosity. We calculate this coupling and show 
that it can in principle significantly modify the stability diagram at lower 
temperatures. As a result, colder stars can remain stable at higher spin 
rates M. 
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In the second one, we propose a possible solution of the unsolved post- 
glitch relaxation of Crab by r-modes. More than 30 years after the discovery 
of the pulsar phenomenon and its identification with neutron stars, there 
exists still a number of uncertainties and open questions about the theoret- 
ical model for pulsars, mainly due to the extremely dense state of matter 
in neutron stars. After two decades, the glitch phenomenon, a sudden in- 
crease of angular velocity of the order of AQ/Q < 10^^, and the very long 
relaxation times, from months to years, after the glitch, remain as one of 
the great mysteries of pulsars. The observed post-glitch relaxation of the 
Crab pulsar has been unique in that the rotation frequency of the pulsar is 
seen to decrease to values less than its pre-glitch extrapolated values. 

The excitation of r-modes at a glitch and the resulting emission of grav- 
itational waves could, however, account for the required "sink" of angular 
momentum in order to explain the peculiar post-glitch relaxation behaviour 
of the Crab pulsar. We show that excitation of the r-modes at a glitch may 
provide a solution to an unsolved observed effect in post-glitch relaxation of 
the Crab pulsar 0. Assuming that r-modes are excited at a glitch, we show 
that this can conveniently describe post-glitch relaxations of both Crab and 
Vela pulsars for a reasonable initial amplitude of the excitation. We use 
a simple model for the total angular momentum of the star, as in |81|, in 



which r-mode amplitude is independent of the rotational frequency of the 
star. 

In chapter 6, we review r-mode instability and CFS mechanism. Further, 
we calculate the coupled shear viscosity-vorticity correction to the r-mode 
timescale. Finally we discuss the possible role of r-mode in post-glitch 
relaxation of the Crab. 
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Chapter 6 

R-mode instabilities in neutron 
stars 

In this chapter wc introduce the recently r-mode instabiUty in neutron stars. 
Thses modes have been found to play an interesting and important role in 
the evolution of hot young rapidly rotating neutron stars. Gravitational 
radiations tend to drive the r-modes unstable in all rotating stars and spin 
down them. 

In section 6.1 we briefly review the r-mode instability. CFS intability, 
the mechanism that governs the r-mode instability, is discussed in section 
6.1.1. In section 6.1.2 the equilibrium model for a slowly and uniformly 
rotating background is elaborated. For small perturbations, the pulsation 
equations of a rotating star in the slow rotation limit are extracted in section 
6.1.3. In Section 6.1.4 the mode eqautions are solved for interesting £ = m 
r-modes. The stability curve of r-mode is discussed in section 6.1.5. 

In section 6.2 using kinetic theory we calculate the effect of shear viscosity- 
vorticity coupling to the stability curve of r-mode. As an application, the 
possible role of r-mode in post-glitch relaxation of the Crab is discussed in 
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section 6.3. 



6.1 r-mode instability 



6.1.1 CFS instability 

The r-mode instability is a member of the class of gravitational radiation 
driven instabilities called CFS instabilities — named for Chandrasekhar, who 



discovered it in a special case and for Friedman and Schutz, who in- 
vestigated it in detail and found that it is generic to rotating perfect fluids 
The CFS instability allows some oscillation modes of a fluid body to 
be driven rather than damped by radiation reaction, essentially due to a 
disagreement between two frames of reference. 

The mechanism can be explained as follows. In a non-rotating star, 
gravitational waves radiate positive angular momentum from a forward- 
moving mode and negative angular momentum from a backward-moving 
mode, damping both as expected. However, when the star rotates the ra- 
diation still lives in a non-rotating frame. If a mode moves backward in 
the rotating frame but forward in the non-rotating frame, gravitational ra- 
diation still removes positive angular momentum — but since the fluid sees 
the mode as having negative angular momentum, radiation drives the mode 
rather than damps it. 

Mathematically, the criterion for the CFS instability is 

uj{uj + mQ)<0, (6.1) 

with the mode angular frequency u (in an inertial frame) in general a func- 
tion of the azimuthal quantum number m and rotation angular frequency 
of star, Q. For any set of modes of a perfect fluid, there will be some modes 
unstable above some minimum m and Q. However, fluid viscosity generally 
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grows with m and also there is a maximum value of Q (known as the Kepler 
frequency Qk) above which a rotating star flies apart. Therefore the insta- 
bility is astrophysically relevant only if there is some range of frequencies 
and temperatures (viscosity generally depends strongly on temperature) in 
which it survives. 

The r-modes are a set of fluid oscillations with dynamics dominated by 
rotation. They are in some respects similar to the Rossby waves found in the 
Earth's oceans and have been studied by astrophysicists since the 1970s . 



The restoring force is the Coriolis inertial "force" which is perpendicular to 
the velocity. As a consequence, the fluid motion resembles (oscillating) 
circulation patterns. The (Eulerian) velocity perturbation is 

(5v = aQR{r/RYYft,e''^' + 0{Q^), (6.2) 

where a is a dimensionless amplitude (roughly 6v/v) and R is the radius 
of the star. is the magnetic type vector spherical harmonic deflned by 

m 

Yf^ = + l)]-i/VV X (rVF,„(^, 0). (6.3) 

Since (5v is an axial vector, mass-current perturbations are large compared 
to the density perturbations. The Coriolis restoring force guarantees that 
the r-mode frequencies are comparable to the rotation frequency []63| , 



uj + mn = ^^n + 0(n^), i = m. (6.4) 
m+1 ^ ^ ^ ^ 



In mid-1997 that Andersson |6^ noticed that the r-mode frequencies satisfy 
the mode instability criterion ( |6.1| ) for all m and Q, and that Friedman and 
Morsink showed the instability is not an artifact of the assumption of 
discrete modes but exists for generic initial data. In other words, a// rotating 
perfect fluids are subject to the instability. 
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6.1.2 Slow rotation approximation 

To analyze the r-modes of rotating stars, we use the standard expansion 
of the hydrodynamics equations as power series in the angular velocity fl 
of the star. In this section we follow the method presented in |1T2|, RSl, 



and describe how to solve the equilibrium structure equations for uniformly 
rotating Newtonian and barotropic stars for slow rotations. The solutions 
will be obtained here up to the terms of order Q^. Here, we use the standard 
spherical coordinates. 

The general equations which describe the dynamical evolution of an ar- 
bitrary state of a Newtonian self-gravitating perfect fluid are the continuity 
equation 

^ + Vaipvn=0, (6.5) 



the Euler's equation 



^ + v'Vhipv") = -V^Mp) - $] = -V"f/, (6.6) 

and the gravitational potential equation 

V"Va$ = -AirGp. (6.7) 

The quantities p, p are the mass density and pressure of the fluid, respec- 
tively. They are assumed to satisfy a barotropic equation of state, p = p{p); 
v"', $, and G are the fluid velocity, the gravitational potential and the 
Newtonian gravitational constant, respectively. Here h{p) denotes the ther- 
modynamic enthalpy of the barotropic fluid in a comoving frame. 

This definition can always be inverted to determine p{h). 
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In equilibrium, we consider a rotating self-gravitating perfect fluid with 
uniform angular velocity, Q. The velocity of the fluid becomes 



v" = (6.9) 

where 0" is the rotational Killing vector field. The equilibrium equations 
will be 

vj/i - _ $] = 0, (6.10) 

V^Va* = -AnGp. (6.11) 

We seek solutions to Eqs. ( |6.10| ) and ( |6.11|) as power series in the angular 
velocity fl. For a slowly rotating star, 

h = ho{r) + 0{n'^), (6.12a) 

p = Po{r) + 0{Q'), (6.12b) 
P = po{r) + O{n^), (6.12c) 
$ = $o(r) + 0(^2)^ (6.i2d) 

where ho, po, po and $0 values for the corresponding non-rotating (spheri- 
cal) equilibrium model. Using these expressions, the zero order solution to 
Eq. (^) is 

Co = ho{r)-%{r), (6.13) 

where Cq is constant. The non-rotating model can be determined in the 
usual way by solving the gravitational potential equation, 

ar \ ar 



together with Eq. ( |6.13|) . The integration constant can be shown to be 



Co = ~GMq/Rq by evaluating Eq. ( |6.13|) at the surface of the star, where 



the constants Mq and Rq are the mass and radius of the non-rotating star. 
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6.1.3 Pulsation equations 

In the last section, the equihbrium model of uniformly rotating star in the 
slowly rotating approximation was discussed. In this section we assume a 
small perturbation in the fluid, to extract the mode equations of the system. 
Using Ipser and Lindblom's method one finds the pulsation equations 
in general. In the next section we restrict our calculations to r-mode only. 

The modes of any rotating barotropic stellar model can be described 
completely in terms of two scalar potentials 5$ and 



where 5p and 5$ are the Eulerian pressure and gravitational potential per- 
turbations, respectively, and p is the unperturbed density of the equilibrium 
stellar model For a barotropic equation of state, p = p{p), Eq. ( |6.15|) 
reduces to 



where Sp is the Eulerian mass density perturbation. We assume here that 
the time dependence of the mode is e*"^* and that its azimuthal angular 
dependence is e**""^, where u is the frequency of the mode and m is an 
integer. 

Linearizing Euler's equation, Eq. ( |6.6| ), about a uniformly rotating back- 
ground, we find 



where Qab is the Euclidean metric tensor (the identity matrix in Cartesian 
coordinates). The quantity can be inverted to obtain an expression for 
the velocity perturbation in terms of the potential 6U : 




(6.15) 



SU = -^6p-6<l> 
p dp 



(6.16) 



iQ-b^6v'' = [i{LU + mn)gab + 2VbVa] = -^JU , 



(6.17) 



(6.18) 
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Q"^^ depends on the frequency of the mode, and the angular velocity of the 
equilibrium star: 



402 

(00 + m9)g''^ -z'^z^ - 2iV%^ 



(6.19) 



In Eq. (|6.19| ) the unit vector 2;" points along the rotation axis of the equi- 
librium star. For real frequencies cj, Q"'^ is Hermitian, = Q* and 
covariantly constant, VcQ"'' = 0. 

Replacing the perturbed mass density and fluid velocity in terms of the 
potentials 6U and 5$, 6p = p{dp/dp){SU + 5$) and using Eq. (|6.18|) , Eqs. 
rsf) and flO) reduce to 



^aipQ^^'^^bSU) = -{iu + mn)p^{SU + 5$), (6.20) 

V«Va<5$ = -A7rGp^{5U + 5$). (6.21) 
dp 

We note that in obtaining Eqs. (|6.20|) and (|6.21|) , the slow rotation ap- 
proximation is not assumed. Equations ( |6.2(]| ) and (|6.21|) are the master 
equations that determine the properties of the oscillations of rapidly rotat- 
ing Newtonian stellar model with uniform spin rate. They are a forth-order 
system of partial differential equations for two potentials 6U and 5$. These 
equations form an eigenvalue problem with eigenfrequencies, u, with appro- 
priate boundary conditions at the surface of the star for 6U and at infinity 
for (5$. 

In the slow rotation limit we expand the potentials 6U and 5$ as 

SU = Rln^6Uo + C»(n^) , (6.22) 



= Rln^S^ + Oin'^). (6.23) 

The normalizations of 6U and 5$ have been chosen to make SUq and 5$o 
dimensionless under the assumption that the lowest order terms scale as 
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Q^. Using Eqs.(|6l9|), (1021) and (lOal) , Eqs. (lOOl) and in the 

lowest order in fl reduce to 



V"U(«o'^"*-4zV)V6^[/o 



'"'°-rV.po 5t/o = 0, 



(5f/o + 5$o). 



(6.24) 



(6.25) 



The quantity Kq is the first order of k, kQ = uj + rriQ. The notation ^ is the 
cylindrical radial coordinate, ^ = rsin('(9), and denotes the unit vector in 
the ^ direction. 

The eigensolutions and eigenfunctions will be determined by solving sys- 
tem of equations ( |6.24| ) and ( |6.25|) with the appropriate boundary conditions 



6.1.4 i = m r-modes 



In this section we restrict our consideration on the r-modes which contribute 
primarily to the gravitational radiation driven instability. The reason for 
restriction to i = m goes back to Provost et al. |]8^, who showed that for the 
barotropic equation of state, only the i = m r-mode exists in the rotating 
star. The "classical" r-modes (which were studied first by Papaloizou and 
Pringle p3]) are generated by hydrodynamic potentials of the form (see e.g. 



Lindblom and Ipser 



6U = a 



^)'^V,+i,(cos(^))e*™'^. (6.26) 



2£ + 1 V ^ + 1 

Here, and after, we drop index for brevity. It is straightforward to verify 
that this b\J is a solution to Eq. ( |6.24| ) for the eigenvalue k has the value 



+ 1' 



(6.27) 
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Then, the frequency of the mode, a;, can be obtained 



(^-l)(^ + 2) 
^ + 1 



n. 



(6.28) 



The perturbed gravitational potential 5$ must have the same angular de- 
pendence as b\J , so 



5$ = a5^(r)Ff+i^(cos(^9))e^'"'^. (6.29) 

Therefore, the gravitational potential Eq. ( |6.25| ) reduces to an ordinary 
differential equation for 5\E'(r): 

dp (£+l)(^ + 2)- 



d^b"^ 2 
+ — 



dr"^ r dr 



AnGp 



dp 



5$ 



dp / r \ 



8nGp£ / 

'2i+l\J+ 1 d \RJ 

(6.30) 



Substituting Eqs. (|6.26| ) and (|6.30|) into Eq. (|6.16| ), the perturbed mass 
density to order fi^ becomes 



p dp 



2i 



2i + iy i+l \R 



(6.31) 



Furthermore, using Eq. ( |6.18|) , the velocity perturbation expression, Eq. 
).2|), will be recovered. 



6.1.5 r-mode timescale 

Our main interest here is to study the evolution of the r-modes due to the 
dissipative influences of viscosity and gravitational radiation. To achieve it, 
we study the energy evolution of the mode which is affected by radiation 
and viscosity. The energy of the mode measured in the rotating frame of 
the equilibrium star, E, is 



E 



p6\ ■ 6\ '' 



' 6p 



(5$ 5p* 



d-^x. 



(6.32) 
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This energy evolves on the secular timescale of the dissipative processes. 
The dissipation of energy due to the gravitational radiation and viscosity 
can be estimated from general expression [p9| 



J 771 



dt 



+ m^]) ^ N,u^' {\5D,„,\'' + |(5J,^p) , (6.33) 

l>2 



where thermodynamic functions 7] and are the shear and bulk viscosities 
of the fluid, respectively. The viscous forces are driven by the shear 5aab 
and the expansion 5a of the perturbation, defined by the usual expressions 

Saab = ^{Va6Vb + VbSVa ' iSabVJv'), (6.34) 

Sa = Va5v\ (6.35) 

Gravitational radiation couples to the evolution of the mode through the 
mass 5D^rn and current 5Jfm multipole moments of the perturbed fluid, 

5Dim = J Sp/Y;^d'x, (6.36) 

SJim = l^J^ J ^'(P^v + 5pv) ■ Yf:^d\ (6.37) 
with coupling constant 



^ 47rG (£+!)(£ + 2) 



The terms in the expression for dE/dt due to viscosity and the gravita- 
tional radiation generated by the mass multipoles are well known |9^. The 



terms involving the current multipole moments have been deduced from the 
general expressions given by Thorne [plj] . 



We can now use Eqs. (|6.32| ) and ( |6.33|) to evaluate the stability of the 



m r-modes. Viscosity, however, tends to decrease the energy while 
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gravitational radiation may either increase or decrease it. The sum that 
appears in Eq. ( 6.33|) is positive definite; thus the effect of gravitational 



radiation is determined by the sign of uioj + IVt). Using Eq. (|6.4|), this 
quantity for r-modes is negative definite: 

u:{u + m) = -^ii^^ll|±^f]2 < 0. (6.38) 

Therefore gravitational radiation tends to increase the energy of these modes. 
In addition, for small angular velocities the energy E is positive definite: 
the positive term |5vp in Eq. ( |6.32| ) (proportional to Q"^) dominates the 
indefinite term {Sp/p — 6^)6p* (proportional to Q"^). Thus, gravitational 
radiation tends to make every r-mode unstable in slowly rotating stars. 
To determine whether these modes are actually stable or unstable in ro- 
tating neutron stars, therefore, we must evaluate the magnitudes of all the 
dissipative terms in Eq. ( |6.33| ) and determine the dominant one. 



Here we estimate the relative importance of these dissipative effects in 
the small angular velocity limit using the lowest order expressions for the r- 
mode 5v and 6p given in Eqs. ( |6.2| ) and ( |6.31| ). The lowest order expression 
for the energy of the mode E is 

E = la^n'^R-^'^^ r P r^'^^dr. (6.39) 
Jo 

The lowest order contribution to the gravitational radiation terms in the 
energy dissipation comes entirely from the current multipole moment 6Ji£. 
This term can be evaluated to lowest order in Q using Eqs. (|6.2|) and (|6.37| ): 



^J» = ^\T^rPr^'^'dr. (6.40) 



ci?^-i V ^ + 1 ^ 

The other contributions from gravitational radiation to the dissipation rate 
are all higher order in Vt. The mass multipole moment contributions are 
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higher order because the density perturbation 5p from Eq. ( |6.31|) is pro- 
portional to Vt^ while the velocity perturbation 5v is proportional to Vt. 
Furthermore, the density perturbation 5p generates gravitational radiation 
at order 2i + A in uj while 5v generates radiation at order 2i + 2. 

The contribution of gravitational radiation to the imaginary part of the 
frequency of the mode I/tg^ can be computed as follows, 



- • (6.41) 



1 _ 1 /dE\ 
Tgr~ 2E\dt)^ 
Using Eqs. (1091) - (|6^) an explicit expression for the gravitational radia- 



tion timescale associated with the r-modes can be obtained: 



The time derivative of the energy due to viscous dissipation is given 
by the shear 6crab and the expansion 6a of the velocity perturbation. The 
shear can be evaluated using Eqs. (|6^) and ( |6.34] ) and its integral over the 



spherical coordinates i) and (p. Using the formulae for the viscous dissipation 
rate Eq. ( |6.33D and the energy Eq. (|6.39|) , the contribution of shear viscosity 



to the imaginary part of the frequency of the mode is, 

— = (£ - 1) (2£ + 1) rjr'^^dr ( p r'^^+^dr] . (6.42) 
Tsv Jo yo J 

The bulk viscosity dissipation expression, 5a, can be re-expressed in 
terms of the density perturbation. The perturbed continuity equation gives 
the relationship 

5a = -i{uj + mfi) Ap/p, (6.43) 

where Ap is the Lagrangian perturbation in the density. The perturbation 
analysis used here is not of sufficiently high order (in Q) to evaluate the 
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lowest order contribution to Ap. However, we are able to evaluate the 
Eulerian perturbation 5p as given in Eq. ( |6.31 ). We expect that the integral 



of |5p/pp over the interior of the star will be similar to (i.e., within about a 
factor of two of) the integral of |Ap/pp. Thus, the magnitude of the bulk 
viscosity contribution to the energy dissipation can be estimated by 

1 {ijj + mVtY f 6p6p 



IC-^d^x. (6.44) 
J 0^ 



Tbv 2E J p' 

Using Eqs. ( |6.31| ) and ( |6.39| ) for Sp/p and E, Eq. ( |6.44| ) becomes an explicit 



formula for the contribution to the imaginary part of the frequency due to 
bulk viscosity. 

To evaluate the dissipative timescales associated with the r-modes using 
the formulae in Eqs. ( |6.42| )-( |6.44| ), we need models for the structures of 



neutron stars as well as expressions for the viscosities of neutron star matter. 
The r-modes timescales for IAMq neutron star models based on several 



realistic equations of state are evaluated by Lindblom et al. [^. The 



standard formulae for the shear and bulk viscosities of hot neutron star 



matter are BS 



7] = 347p^/^T-^gcm-^s-\ (6.45) 

C = 6.0 X 10~^V^(^ + mn)-^T^gcm-h-\ (6.46) 

The timescales for the more realistic equations of state are comparable 
to those based on a simple polytropic model p = kp^ with k chosen so that 
the radius of a IAMq star is 12.53 km. The dissipation timescales for this 
polytropic model (which can be evaluated analytically) are Tq^ = ~ 3.26s, 
fsv = 2.52 X lO^s and fgy = 6.99 x lO^s for the fiducial values of the angular 
velocity Q = ^pnUp and temperature T = lO^K in the £ = 2 r-mode. Here 
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p = 3M/AttR^ is the mean density of the star. The gravitational radiation 
timescales increase by about one order of magnitude for each incremental 
increase in i, while the viscous timescales decrease by about 20%. 

The evolution of an r-mode due to the dissipative effects of viscosity and 
gravitational radiation reaction is determined by the imaginary part of the 
frequency of the mode, 



' ^ +^ T7^ -7T-J- (6-47) 



t{^) Tgr xT^Gp) Tsv \ T J Tbv VlO^K/ \7tGp^ 

Eq. (|6.47 ) is displayed in a form that makes explicit the angular velocity 
and temperature dependences of the various terms. Dissipative effects cause 
the mode to decay exponentially as e"*/"^ (i.e., the mode is stable) as long 
as r > 0. From Eqs. ( |6.42| )-( py4^ ) we see that fsv > and > while 
T^GR < 0. Thus gravitational radiation drives these modes towards insta- 
bility while viscosity tries to stabilize them. For small Q the gravitational 
radiation contribution to the imaginary part of the frequency is very small 
since it is proportional to fi^'"*"^. Thus for sufficiently small angular veloc- 
ities, viscosity dominates and the mode is stable. For sufficiently large Q, 
however, gravitational radiation will dominate and drive the mode unsta- 
ble. It is convenient to define a critical angular velocity Qc where the sign 
of the imaginary part of the frequency changes from positive to negative: 
1 / t{Qc) = 0. If the angular velocity of the star exceeds Qc then gravitational 
radiation reaction dominates viscosity and the mode is unstable. 

For a given temperature and mode / the equation for the critical angular 
velocity, = 1/t{Qc), is a polynomial of order / + 1 in Q^, and thus each 
mode has its own critical angular velocity. However, only the smallest of 
these (always the / = 2 r-mode here) represents the critical angular velocity 



of the star. Fig. 6.1 shows that the critical angular velocity for a range of 
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temperatures relevant for neutron stars for the polytropic model discussed 



above. Fig. 3.2 depicts the critical angular velocities for IAMq neutron star 



models computed from a variety of realistic equations of state ||92[ . Fig. 



illustrates that the minimum critical angular velocity (in units of ^JirGp) is 
extremely insensitive to the equation of state. The minima of these curves 
occur at T ^ 2 X lO^K, with VLc ~ 0.043A/7rG'p. The maximum angular 
velocity for any star occurs when the material at the surface effectively 



orbits the star. This 'Keplerian' angular velocity Vt^ is very nearly ^y/nGp 
for any equation of state. Thus the minimum critical angular velocity due 
to instability of the r-modes is about 0.065flK for any equation of state. 
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Figure 6.1: Critical angular velocities for a IAMq polytropic neutron star. 
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Figure 6.2: Critical angular velocities of realistic IAMq neutron star mod- 
els, Lindblom et al. IBSl 
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6.2 Vorticity-shear viscosity coupling 



In spite of the recent improvements in our understanding of the r-mode 
instabihty, it seems that the fundamental properties of these modes have 
not yet been sufficiently understood. Previous investigations of the r-modes 
are restricted to the case of uniformly and slowly rotating, isentropic, New- 



tonian stars ||S^. A few recent studies were done for relativistic stars with 
slowly rotating and Cowling approximations In this sense, it is inter- 
esting to study the properties of the r-mode instability in the more general 
cases, for example, differentially and rapidly rotating, non-isentropic rela- 
tivistic stars. 

In addition, people have used the standard Navier-Stokes theory to study 
r-mode stability enforced by viscosity, and have calculated the correspond- 
ing timescales. This theory and its relativistic generalization are non-causal 
and unstable. An improved causal dissipative fluid theory is based on kinetic 
theory [|l^. This latter theory has implications different from Navier-Stokes 
theory. For example, kinetic theory predicts that the angular velocity of the 
star couples to the viscosity and heat flux in rotating stars. 

In this section we investigate the possible effect of vorticity-shear vis- 
cosity coupling on the stability of r-mode and show that the coupling may 
have a significant effect on viscous damping timescales of r-mode Q]. We 
find that colder stars can remain stable at higher spin rates. 

6.2.1 Predicted effect 

In standard Navier-Stokes theory the viscous quantities are defined by 

n = -C© , ga = -/tV^T , -Kab = -2vcrab , (6.48) 
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where 11 is the bulk viscous stress and ( is the bulk viscosity; qa is the 
heat flux and k is the thermal conductivity; TTab is the shear viscous stress 
and 7] is the shear viscosity; = VaV"" is the volume expansion rate of the 
fluid, T is the temperature, and aab = ^{aVt) is the rate of shear (where 
the angled brackets denote the symmetric tracefree part). Its is clear that 
the fluid vorticity tu = |V x v does not enter Eq. ( |6.48D , even when the 
equilibrium state is rotating. 

On physical grounds, one might expect that rotational accelerations can 
couple with gradients of momentum and temperature, so that there could 
in principle be couplings of Wa to qa and iTab- In the case of heat flux, 
qualitative particle dynamics indicates ||9^ (p. 34) that this coupling does 
exist as a result of a Coriolis effect, which is in some sense analogous to the 
Hall effect in a conductor subject to a magnetic field. The Coriolis effect on 
heat flux is confirmed by molecular dynamics simulations |Q. Miiller 



and Israel & Stewart |]T8[ showed that the Boltzmann equation predicts in 
general a coupling of vorticity to heat flux and shear viscous stress. The 
microscopic and self-consistent kinetic approach is in contrast to the con- 
tinuum view, where a phenomenological principle of "frame indifference" is 
invoked to argue against any vorticity coupling. (See ^ for further 

discussion.) 

Using the Grad moment method to approximate the hydrodynamic 
regime via kinetic theory, the relations in Eq. ( |6.48|) are modified to [jl8| 
(Eq. (7.1)) 



n 

Qa 
T^ab 



-C e + /?on 



-K 

-2r] 



VaT + m [qa - -^abq'} 



(6.49a) 
(6.49b) 
(6.49c) 
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where (3 a-, A = 0,1,2, can be evaluated in terms of collision integrals for 
specific gases, an overdot denotes the comoving (Lagrangian) derivative, 
and the vorticity tensor is given by 



w = X V , 

where square brackets on indices indicate the skew part. Navier-Stokes 
theory is recovered from the Miiller-Israel- Stewart theory when (3a = 0. 
However, kinetic theory gives (3a values for simple gases which are definitely 
not zero. Furthermore, if (3a = 0, the equilibrium states are unstable and 
dissipative signals can propagate at unbounded speed |]18, 



The /3^-corrections will be very small except if there are either high fre- 
quency oscillations (pumping up the time- derivative terms) or rapid rotation 
(pumping up the vorticity-coupling terms). In the context of rapidly rotat- 
ing neutron stars, we expect the vorticity-dissipative couplings to dominate 
the time-derivative terms; this expectation is borne out by calculations (see 
below). The vorticity-dissipative couplings will be negligible if the unper- 
turbed equilibrium state is irrotational, i.e., if zua = in the background, 
so that the coupling terms become second-order. However, for fast rota- 
tion, Wa 7^ in the background and the coupling terms make a first-order 
contribution to dissipation. In the words of Israel & Stewart |18[: "these 
results will ultimately be of practical interest in astrophysical and cosmo- 
logical situations involving fast rotation, strong gravitational fields or rapid 
fluctuations (neutron stars, black hole accretion, early universe), although 
it will probably be some time before the state of the art in these fields makes 
such refinements necessary." 

We believe that recent and ongoing developments in rotating neutron star 
physics have reached the stage where the Miiller-Israel- Stewart theoretical 
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corrections to the Navier-Stokes equations need to be examined, and our 
results indicate that the corrections could be important. 

We follow the standard assumption that the heat flux may be ne- 
glected relative to viscous stresses in calculating damping timescales. Then 
the vorticity correction to Navier-Stokes theory reduces to the coupling term 
^^(a7i"fe)c- This term means that the angular momentum of the star changes 
the shear viscosity timescale, and we find (for axial r-modes) a correction 
proportional to T^^Q"^, where r = 9 for a nonrelativistic fluid and r = 12 
for an ultrarelativistic fluid. 

Replacing 6(7 and 6aab by 611 and 6Trab in Eq. ( |6.33D , respectively, the 
evolution of dissipation energy contained in small fluctuations is given by 

dE f f6mu* 67T''''6TT*f,\ 3 fdE 



,d^x-\ — \, (6.50) 
dt J \ C 2r^ J \dt J^^^ ^ ' 

where {dE / dt)^^ is the energy flux in gravitational radiation (see Eq. ( 6.33| )), 



5n = n — n and 6Tiab = T^ab — T^ab, with an overbar denoting background 
quantities. In this case, n = = Ttab- The normal modes of the star are 
damped by dissipation, and the damping rate can be determined by Eq. 
(|6.50| ). For a normal mode with time dependence e*"^*, the energy has time 
dependence exp[—2Im(ti;)t]. Then by Eq. ( |6.50|) , the characteristic damping 
time r = l/lm{uj) of the fluid perturbation is given by 

1 1 dE 1 1 1 

- = --77 = — + — + — ' (6-51) 

T 2E dt Tbv Tgv Tqr 

where Tgy, Tgv, and Tqr are the bulk viscous, shear viscous, and gravitational 
radiation timescales respectively. 

To evaluate the vorticity-corrected shear viscous timescale, we use Eq. 
( |6.49c|) in Eqs. ( |6.50|) and ( |6.51|) . To lowest order 



S'n'ab = -2?7 6aab - 2iujr](326aab + 4r](326(j?^Wh) 
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where Wa is the background vorticity (the background shear vanishes) . Then 



where 7 = ?7/32- The first term is the usual term in Navier-Stokes theory, 
while the following terms are the Miiller-Israel- Stewart corrections. The u'^ 
term arises from -kab in Eq- ( |6.49c| ), and is negligible relative to the w"^ terms 
which arise from the w'^/^a'^^'jc term in Eq. ( |6.49c|) . The energy dissipation 
rate through shear viscosity will be 



+ 4 ((5(T"^5a*bt<7"t<7e - } d^x. (6.52) 

In order to proceed further, we need expressionx for the shear viscosity 
7] and the coupling coefficient /52- For the various interactions, ri{p,T) is 
calculated in [0, where it is shown that electron-electron scattering is 
more important for shear viscosity than other interactions. The expression 
for 7] is given in [R^, in good agreement with p^, p8 



as 



7] = 1.10 X 10 



16 



T 



g/cms . 



.53) 



For a Maxwell-Boltzmann gas, the coefficient (32 is found in but we 
require the expression for a degenerate Fermi gas. This has been found by 
Olson & Hiscock in the case of strong degeneracy: 



P2 



15nm-' (1 + z/) 
m^gc^ (z/2 + 2z/)5/2 



+ 



' kT ' 
mc^u 



(6.54) 
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where m is the particle mass, g is the spin weight, and mc^v/kT ^ 1. 
The dimensionless thermodynamic potential u = {p + p)/nm — mc^s/kT — 
1, where s is the specific entropy, is equal to the nonrelativistic chemical 
potential per particle divided by the particle rest energy. For a strongly 
degenerate gas, the nonrelativistic chemical potential is proportional to T, 
so that 

_ kT 
u ^ a 



where a ^ 1 is a dimensionless constant measuring the degree of degener- 
acy. The nonrelativistic regime is obtained for z/ <^ 1, while the ultrarela- 
tivistic case corresponds to z/ ^ 1. 

For temperatures below 10^° K, neutrons in the neutron star are nonrel- 
ativistic, while electrons are ultrarelativistic The nonrelativistic limit 
of P2 is 

{P2)nr ~ 3.16 X 10-\aT)-^/^ cm sVg , (6.55) 
and its ultrarelativistic limit is 

(/?2)uR, ~ 6.45 X lO^^aT)-^ cm sVg • (6.56) 
Using Eqs. (|6.53|) , ( |6.55|) and ( |6.56D , we have 



1.10 X 10-^v p y /Vio^Ky/- 



7.08 X 10-5 / p flO^K\ 



(6.58) 



VlQi^g/cmV \ T J ' 

In these calculations, we have used the same relation for rj in both cases, 
because in the high-density regime (p > lO^^g/ cm^) for both electron- 
electron scattering and electron-neutron scattering, rj is proportional to 
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T ^, with nearly equal proportionality factor For typical values of 

the temxerature, T = 10^ K, and density, p = 3 x 10^*^ g/cm^, we find that 
7uR ~ a-^ X 10"^ s, while ~ a'^'^ x 10-^° s. 

6.2.2 r-mode instability curve 



In this section we calculate the predicted effect, Eq. ( 6.52| ), for the r-mode 



instability. We assume that the background is a uniformly rotating star, so 
that the equilibrium fluid velocity is f° = VLcf)'^, where cj)"" is the rotational 
Killing vector field ||10CI|| . The vorticity vector of the equilibrium state is 

w=^ [cot ^9, -1,0] . (6.59) 

The r-modes of rotating barotropic Newtonian stars have Eulerian velocity 
perturbations given by Eq. ( |6.1^ ) 

5v = am Yff, exp{iut) , (6.60) 

where C is an arbitrary constant, R is the unperturbed stellar radius, and 
u = 2^2/(^+1). The magnetic-type vector spherical harmonics are 
defined by 



^Fm = , , V X [rVr,™(t9, ^)] . (6.61) 
£+1 



The shear of the perturbed star is given by 



5(Tab = V(a(5w6) ■ (6.62) 



Substituting Eqs. ( |6.5!j| )-( |B.62| ) into Eq. ( p.52| ), we find the shear viscosity 



timescale for i = m: 
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(i - 1) (2£ + 1) r]r^^ dr + il^S, 
Jo 



(6.63) 



where Qi^ = pr'^^'^'^ dr . The first term in brackets is in agreement with 



the expression calculated in |77|, and Si is the correction term: 



16 



(£-i)(2£ + i)^^ , i{i - 2)\[{2i - ly.if r(i) 



{£ + iy ' {i + i){2i - i){2i)\ r(£-i; 

{2f - 8f -U- 6)f/2 + I2{e -f-e+ l)f/3 



+ 2{4i^ - f - 9e^ + 5e + i)U4 

where f/fc(T) = R'^J^-f^ 7]r^^-^dr. 

For the ^ = 2 modes, Eqs. ( |6.63|) and (|6.64|) give 



(6.64) 



1 = 5Q2 r dr + \Q2^'' [80?7o + 93[/2 + bW^ - A2Ui] . (6.65) 

Ts Jo 

For comparison with previous calculations based on Navier-Stokes viscosity 
(see, e.g., |8^), we use an n = 1 polytrope with mass M = IAMq and radius 
R = 12.53 km to evaluate the integrals in Eq. ( |6.65| ). The bulk viscous and 
gravitational radiation timescales are unaffected by the vorticity correction, 
and we obtain 

1 1 r 1 / T \ V \^ 



r(fi, T) Va^J ^ \Wk 

Tsv V T I 



qa 



4-r 



T 



(6.66) 



where f^K = y/irGp, which is | times the Keplerian (mass-shedding) fre- 
quency, and the vorticity correction factors are 



1.36 X 10 
5.67 X 10 



-23 



-10 



9 nonrel , 
12 ultrarel . 



(6.67) 
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The standard result (see, e.g., [^) is regained for q = 0, with 



f^^ = -3.26 s , fev = 2.01 X 10^^ s , fsv = 2.52 x 10" s . 

We note that the contribution from the -kab term in Eq. ( |6.49c| ) to the 
g-correction is less than 1% of the contribution from the cu^^aVr^^c term. 



Now we are able to determine from Eq. (|6.66|) the critical angular ve- 
locity defined by 1/t{Qc,T) = 0, which governs stability of the star: 
if > Qc, then dissipative damping cannot overcome the gravitational 
radiation-driven instability. In Fig. |6.3| we plot Qc/^k against temperature 
T, showing how the vorticity-viscosity coupling affects the standard result 
(see, e.g., [1^). Electrons are assumed to dominate the shear viscosity, and 
they are ultrarelativistic over the range of temperatures. 



It is clear from Fig. ^]3| that the vorticity correction is only appreciable 
at temperatures T < 10^ K, but that for these lower temperatures, the 
correction can be large, especially for smaller a. As the degree of degen- 
eracy increases (i.e., with increasing a), the correction is confined to lower 
and lower temperatures. The effect of the vorticity-viscosity coupling is to 
increase the stable region, so that cooler stars can spin at higher rates and 



remain stable. This may modify recent results [[L01|] which suggest that 
r-mode instability could stall the spin-up of accreting neutron stars with 
T > 2 X 10^ K; if the vorticity correction operates, then the stability region 
is increased, so that spin-up could be more effective, especially for lower 
degeneracy parameter a. 

We note that, here, we have used for our r-mode calculations solutions 
that assume slow rotation. Thus the Q/Qj^ > 0.3 part of Fig. |6.3| is an ex- 
trapolation to high spin rates, in common with previous stability diagrams. 
Recent calculations of r-modes for rapid rotation should be used in fu- 
ture calculations of the vorticity correction. Since /-modes are unstable at 
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high spin rate, the effect of the vorticity correction on these modes would 
also be interesting to calculate. 



85 




T(xlO^) K 

Figure 6.3: Critical angular velocity versus temperature {n — 1 polytrope 
with mass IAMq and radius 12.57 km). The stability region is below the 
curves. The solid curve shows the standard result, with no coupling of 
viscosity to vorticity. Broken curves (labelled by the degeneracy parameter 
a) show how the instability region is reduced by the kinetic-theory coupling 
of shear viscosity to vorticity, for an ultra-relativistic degenerate Fermi fluid 
(electron-electron viscosity). 



86 



6.3 Post-glitch relaxation of Crab 

As we discussed before, the r-mode instability opens a wide window for 
gravitational wave astronomy. It would give us some information about the 
interior matter of neutron stars. Determining cooling rates, viscosity, crust 
formation, the equation of state of neutron matter, the onset of superflu- 
idity in neutron stars, and several other features of neutron stars are more 
interesting implications of this instability. In this section, we discuss one of 
possible r-mode implications. 

More than 30 years after the discovery of the pulsar phenomenon and its 
identification with neutron stars, there exists still a number of uncertain- 
ties and open questions about the theoretical model for pulsars, mainly due 
to the extremely dense state of matter in neutron stars. During the past 
two decades, the glitch phenomenon, a sudden increase of angular velocity 
of the order of An/n < 10'^, and the very long relaxation times, from 
months to years, after the glitch, remain as one of the great mysteries of 
pulsars. The observed post-glitch relaxation of the Crab pulsar has been 
unique in that the rotation frequency of the pulsar is seen to decrease to 
values less than its pre-glitch extrapolated values. So far, two mechanisms 
have been suggested to account for the observed excess loss of angular mo- 
mentum during post-glitch relaxations of the Crab. The first mechanism, 



in the context of the vortex creep theory of Alpar et al. | 102 , 103], in- 
vokes generation, at a glitch, of a so called "capacitor" region within the 
pinned superfiuid in the crust of a neutron star, resulting in a permanent 
decoupling of that part of the superfiuid. Neverthless, this suggestion has 
been disqualified (Lyne et al. ||104|| ) since the moment of inertia required to 



have been decoupled permanently in such regions during the past history 
of the pulsar is found to be much more than that permitted for all of the 
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superfluid component in the crust of a neutron star. In another attemp, 



Link et al. ||105|| have attributed the excess loss of angular momentum to 
an increase in the electromagnetic braking torque of the star, as a conse- 
quence of a sudden increase, at the glitch, in the angle between its magnetic 
and rotation axes. As they point out, such an explanation is left to future 
observational verification since it should also accompany other observable 
changes in the pulsar emission, which have not been detected, so far, in 
any of the resolved glitches in various pulsars. Moreover, the suggestion 
may be questioned also on the account of its long-term consequences for 
pulsars, in general. Namely, the inclination angle would be expected to 
show a correlation with the pulsar age, being larger in the older pulsars 
which have undergone more glitches. No such correlation has been deduced 
from the existing observational data. Also, and even more seriously, the 
assumption that the braking torque depends on the inclination angle is in 
sharp contradiction with the common understanding of pulsars spin-down 
process. The currently inferred magnetic field strengths of all radio pulsars 
are in fact based on the opposite assumption, namely that the torque is in- 
dependent of the inclination angle. The well-known theoretical justification 
for this, following Goldreich & Julian | [LO(j|| , is that the torque is caused by 
the combined effects of the magnetic dipole radiation and the emission of 
relativistic particles, which compensate each other for the various angles of 



inclination (see, eg., Manchester & Taylor ||107|| ; Srinivaran ||108 



The excitation of r-modes at a glitch and the resulting emission of grav- 
itational waves could, however, account for the required "sink" of angular 
momentum in order to explain the peculiar post-glitch relaxation behav- 



ior of the Crab pulsar. As is shown in Figs. and |6.5| , for values of 
ctQ > 0.04 the predicted time evolution of ^ and ^ during the 3-5 years 
of the inter-glitch intervals in Crab, might explain the observations. That 
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is, the predicted total change in the rotation frequency of the star, | | , is 
much larger than the corresponding jump ^ ~ 10~^ at the glitch, which 
explains why the post-glitch values of Q should fall below that expected 
from an extrapolation of its pre-glitch behavior. Also, the predicted values 
of ^ ~ 10~^, after a year or so (Fig. 6^), are in good agreement with 
the observed persistent shift in the spin-down rate of the Crab (Lyne et al. 



The predicted increase in the spin-down rate would be however dimin- 
ished as the excited modes at a glitch are damped out, leaving a permanent 
negative offset in the spin frequency. Hence the above so-called persistent 
shift in the spin-down rate of the Crab may be explained in terms of the 
effect of r-modes, as long as it persists during the inter-glitch intervals of 2-3 
years. It may be noted that a really persistent shift in the spin-down rate 
at a glitch may be caused by a sudden decrease in the moment of inertia 
of the star. However this effect, by itself, could not result in the observed 
negative offset in the spin frequency. 

The same mechanism would be expected to be operative during the post- 
glitch relaxation in the other colder and slower pulsars, as well. However, 
for the similar values of Oq, ie. the same initial amplitude of the excited 
modes, the effect is not expected to become "visible" in the older pulsars. 
Particularly, for the Vela its initial jump in frequency at a glitch, ^ ~ 
10~^, is seen from Fig. |6.4| to be much larger (ie. by some four orders of 
magnitudes) than that of the above effect due to the r-modes. In other 
words, while the predicted loss in the stellar angular momentum due to the 
excitation of r-modes result in a negative ^ which, in the case of Crab, 
overshoots the initial positive jump at a glitch, however for the Vela and 
older pulsars it comprises only a negligible fraction of the positive glitch- 
induced jump. A more detailed study should, however, take into account 



89 



the added complications due to internal relaxation of various components 
of the star, which is highly model dependent. The observed initial rise in 
Q need not be totally compensated for by the losses due to r-modes which 
we have discussed, since part of it could be relaxed internally (by a transfer 
of angular momentum between the "crust" and other components, and/or 
temporary changes in the effective moment of inertia of the star) even in 
the absence of any real sink for the angular momentum of the star. Such 
considerations would not only leave the above conclusions valid but also 
allow for even smaller values of the initial amplitude of the excited modes, 
compared to our presently adopted value of ao ~ 0.04. 

The suggested effect of the r-modes in the post-glitch relaxation of pul- 
sars should be understood as one operating in addition to that of the internal 
relaxation which is commonly invoked. While the latter could account only 
for a relaxation back to the extrapolated pre-glitch values of the spin fre- 
quency, the additional new effect due to the r-modes may explain the excess 
spin-down observed in the Crab pulsar as well. 

It is further noted that the above estimates are for an adopted value 
of Q = 9.4 X 10^^, which corresponds to the particular choice of the poly- 
tropic model star. Differences in the structure among pulsars, in particular 
between Crab and Vela, which have also been invoked in the past (see, eg., 
Takatsuka & Tamagaki [|110|| ), could be further invoked to find a better 
agreement with the data for the above effect due to r-modes as well. Also, 
the initial amplitude of the excited modes need not be the same in all pul- 
sars. It is reasonable to assume that in a hotter and faster rotating neutron 
star, as for the Crab, larger initial amplitudes, ie. larger values of ao, are 
realized than in the colder-slower ones. 
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6.3.1 Predicted Effect 



In order to estimate the effect of tlie r-mode instability, in a "stable" neu- 
tron star, on its post-glitch relaxation, we have used the model described by 



Owen et al. ||8T|- The total angular momentum of a star is parameterized 
in terms of the two degrees of freedom of the system. One, is the uniformly 
rotating equilibrium state which is represented by its angular velocity Qeq- 
The other, is the excited r-mode that is parameterized by its magnitude a 
which is bound to an upper limiting value of a = 1, in the linear approxi- 
mation regime treated in the model. Thus, the total angular momentum J 
of the star is written as a function of the two parameters fleq and a: 

J = Jeqfieq + Jc, (6.68) 



'eq"^ "eq 

where I^q = IMR^ is the moment of inertia of the equilibrium state, and 



Jc = —^Ja'^QeqMR'^ is the canonical angular momentum of the I = 2 r- 
mode, which is negative in the rotating frame of the equilibrium star. The 
dimensionless constants 



in 



^=^2L P^'dT (6.69) 







1 

J = — — / pr^dr, (6.70) 



depend on the detailed structure of the star, and for the adopted n = 1 
polytropic model considered have values / = 0.261 and J = 0.01635. Also 
R = 12.54 km and M = 1.4 Mq are the assumed radius and mass of the 
star, for the same polytropic model. 

Eq. ( p.68| ) above implies that an assumed instantaneous excitation of 
r-modes at a glitch would cause a sudden increase in Q^q- For definite- 
ness, we define the "real" observable rotation frequency Q of the star as 



91 



= J, where I is the moment of inertia of the real star. The two are equal, 
Q = Qeq, in the absence of the r-modes, ie. before the excitation of the 
modes at a glitch and after the modes are damped out. If there were no 
loss of angular momentum (by gravitational radiation) accompanying the 
post-glitch damping of the modes (by viscosity) Q^q would recover its ex- 
trapolated pre-glitch value; ie. its initial rise would be compensated exactly. 
However due to the net loss of angular momentum by the star, the post- 
glitch decrease of Q^q overshoots its initial rise. The negative offset between 
values of Qeq before the excitation of the modes and after they are damped 
out is the quantity of interest for our discussion. The question of whether 
the instantaneous rise in the value of Q^q at a glitch, due to the excitation 
of the r-modes, is observable or not is a separate problem, and its resolution 
would have no consequence for the net loss of angular momentum from the 
star which is the relevant quantity here. It is noted that the distinction 
between fl and fleq, in the presence of modes, is quantitatively negligible, in 
all cases of interest, and is usually disregarded. Also, one might dismiss an 
increase in as implied by Eq. ( |6.68|) to be observable as a spin-up of the 



star since for an inertial outside observer the r-modes rotate in the prograde 
direction and their excitation should result, if at all, in a spin-down of the 
star. Moreover, an excitation of the r-modes should not result, by itself, in 
any real change of the rotation frequency of the star at all. Because one 
could not distinguish two physically separate parts of the stellar material 
such that the two components of angular momentum in Eq. (|6.68|) may be 
assigned to the two parts separately. 

The total angular momentum J of the star in terms of Qeq and a is 

J(fieq, a)=(i- pa^] Q,qMR\ (6.71) 



2 

The perturbed star loses angular momentum primarily through the emis- 
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sion of gravitational radiation. Thus, the evolution of J(fieq, «) can be com- 
puted by using the standard multipole expression for angular momentum 
loss. Then, for the i = m = 2 case 



(322)', (6.72) 



dt levrG V 3 

where c is the speed of light and i = m = 2 current multipole, ^22, is given 
by [H 



^22 = V2^^an,^R'J. (6.73) 
15 



Combining Eq. ( 6.72| ) for the angular momentum evolution of the star with 



Eqs. ( |6.42| ) and ( |6.73| ), we obtain one equation for the evolution of the 
parameters and a that determine the state of the star: 

In addition to radiating angular momentum from the star via gravita- 
tional radiation, the mode will also lose energy through gravitational radi- 
ation and neutrino emission (from bulk viscosity). Furthermore the mode 
energy is deposited into the thermal state of the star by shear viscosity. 
Therefore the energy balance equation should be considered together with 
Eq. ( |6.74|) to determine the parameters f2eq and a. For the i = 2 r-mode 
E, Eq. (|6.39|) , is given by 

E = la^nl^MR^J. (6.75) 
The time derivative of E, Eq. ( |6.51D , is 

^4 = -2E f 1 + ^1 . (6.76) 



dt \ Ty Tgr 
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where Tv and Tqr are the viscous and gravitational radiation timescales re- 
spectively. Combining Eqs. ( |6.75| ) and (|6.76|) , the second evolution equation 
for f2eq and a is given by 



da 



da 



rf^ + " dt 



eq 



-afieq (- + — ). (6.77) 

Therefore the time evolution of the quantities a and can be deter- 
mineded from the coupled equations (Owen et al. jSlI): 



da a a 1 — o?Q 



dt Tgr Tv 1 + a^Q ' 



(6.78b) 



where Q = \j = 0.094, for the adopted equilibrium model of the star. The 
viscous time has two contributions from the shear and bulk viscousities with 
corresponding timesacels and Tbv, respectively. The overall "damping" 
timescale r for the mode, which is a measure of the period over which the 
excited mode will persist, is defined as 

111111 , ^ 

- = — + — = — + — + — (6.79) 

'^v ''"gr, 'HjV ''"bv ''"gr 

Following Owen et al. we use Tgv = 2.52 x 10^(s)rg^ Tbv = 4.92 x 
10^°(s)T9"^fi3^ and Tgr, = -1.15 x 10^(s)fiJ^ where Tg is the temperature, 
T, in units of 10^ K, and ^3 is in units of lO^rad s~^. These estimates do 
not however include the role of superfluid mutual friction in damping out 
the oscillations. We have further taken into account the damping due to 
the mutual friction using the associated damping time as given by Lindblom 



and Mendell ||111|| . The effect of the mutual friction is nevertheless seen to 
be negligible and the computed curves shown below remain almost the same 
in the presence of mutual friction. 
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By integrating Eqs. ( |6.78a]) and (|6.78b|) , numerically, for a given initial 
value of a, one may therefore follow the time evolution of a and fieq which 
together with Eq. ( |6.68| ) determine the time evolution of the total angular 
momentum, J, and hence the time evolution of fi. Figs. |6.4| and |6.5| shows 
the computed time evolution for the absolute value of the resulting (neg- 
ative) fractional change ^ in the spin frequency (Fig. |6.4D and also the 
change ^ in the spin-down rate of the star (Fig. |6.5|) , starting at the glitch 
epoch which corresponds to time t = 0. The results in Figs. p]4| and B]5| are 



for a choice of an initial value of ao = 0.04, and for the assumed values of 
T and corresponding to the Crab and Vela pulsars, as indicated. Fig. |0 
shows that for the same amplitude of the r-modes assumed to be excited 
at a glitch the resulting loss of angular momentum through gravitational 
radiation would be much larger in Crab than in Vela, ie. by more than 3 
orders of magnitudes. (Note that the curve for Vela in Fig. ^!4| represents 
the results after being multiplied by a factor of 10^.) Furthermore, for the 
adopted choice of parameter values, the magnitude of the corresponding 
decrease in Vl for the Crab, is |^| ~ 10^^ (Fig. |6^) . The observational 
consequence of such an effect would nevertheless be closely similar to what 
has been already observed during the post-glitch relaxations of, on/y, the 
Crab pulsar. 

Before proceeding further with Crab, we note that the post-glitch effects 
of excitation of r-modes would however have not much observational con- 
sequences for the Vela, and even more so for the older pulsars, which are 
colder and rotate more slowly. This has two, not unrelated, reasons: in the 
older pulsars r-modes a) are damped out faster (ie. have smaller values of 
r), and b) result in less gravitational radiation. The dependence of r on the 



stellar interior temperature is shown in Fig. |6.6| . For the colder, i.e. older, 
neutron stars the r-modes are expected to die out very fast. The damping 
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timescale for a pulsar with a period P ~ Is, being colder than 10^ K, could 
be as short as a few hours (Fig. |6.6| ), and r- modes would have been died 
out at times longer than that after a glitch. For the hot Crab pulsar, on 
the other hand, r-modes are expected to persist for 2-3 years after they are 
excited, say, at a glitch. The value of r decreases for older pulsars due to 
both their longer periods as well as lower temperatures, but the effect due 
to the latter dominates by many orders of magnitudes, for the standard 



cooling curves of neutron stars (Urpin et. al. [|112|| ). The second reason, ie. 
the loss of angular momentum being negligible in older pulsars, was already 
demonstrated in Fig. |6.4| , by a comparison between Crab and Vela pulsars. 
We have verified it also for the case of pulsars older than Vela. It may be 
also demonstrated analytically from Eqs. (|6.78a| ) and (|6.78b| ), in the limit 
of a^Q << 1. The initial increase in f2eq due to excitation of r-modes with a 
given initial amplitude ao is seen from Eq. to be |^|o = alQ. The 

subsequent damping of the modes result in secular decrease in f2eq, and the 
total decrease at large t ^ oo would be l^p^loo ~ — «nQ) which is true for 

i icq "^V ^ 

lA^eql << fieq- Notc that in the absence of gravitational radiation losses 
(ie. — = 0; r = Tv) the total decrease would be the same as the initial 

^ Tgb ' ' 

increase, which is expected for the role of viscous damping alone. The dif- 
ference between these two changes (total decrease minus initial increase) in 
f2cq would correspond to the total loss of angular momentum from the star, 
hence to the net decrease in its observable rotation frequency, ie. 



r«o<5, (6-80) 



which is valid in the limit of << 1. Fig. |6.7| shows the dependence of 
the quantity on the stellar rotation frequency, and also on its internal 
temperature. While for the Crab ~ 10~^, however its value is much 
less for the older pulsars, due to both their lower Vl as well as lower T 
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values. The dependence on the temperature is however seen to be much 
less than that on the rotation frequency, in contrast to the dominant role 
of the temperature in determining the value of the total damping time r, 
as indicated above. As is seen in Fig. 16. 71, for the Vela < 10~^, which 
means a maximum predicted value of |^|oo < 10~^, even for the large 
values of ao ~ 1. This has to be contrasted with the glitch induced values 
of ~ 10^® in Vela, which shows the insignificance of the role of r-modes 
in its post-glitch behaviour. 
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Figure 6.4: The post-glitch time evolution of the spin frequency of a pulsar, 
caused by its loss of angular momentum due to gravitational waves driven by 
the r-modes that are assumed to be excited at the glitch epoch, t = 0, with 
an initial amplitude of = 0.04. The two curves correspond to assumed 
values of Tg = 0.3 and ils — 0.19, for the Crab {thickline), and Tg — 0.2 and 
— 0.07, for the Vela {thin hne). Note that the curve for Vela represents 
the results after being multiplied by a factor of 10^. 
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Figure 6.5: Time evolution of the spin-down rate of a pulsar, caused by its 
loss of angular momentum due to the excitation of r-modes at t = 0. A 
value of 17 = 2.4 x 10~^rad s~^, and other parameter values same as in Fig. 
|6.4| for the Crab have been assumed. 
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Figure 6.6: Dependence of the total damping timescale of r-modes on the 
internal temperature of a neutron star. Parameter values same as for the 
Crab in Fig. |6.4| have been assumed. 
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Figure 6.7: Dependence of the net post-glitch decrease in the rotation 
frequency, on the rotation frequency of a neutron star. The two curves are 
to show the dependence on the temperature, where Tg = 0.3 (bare hne) and 
Tg = 0.2 {dotted hne) have been used. 
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Chapter 7 



Concluding remarks 



In many cosmological and astrophysical situations, an idealized fluid model 
of matter is inappropriate, and a self-consistent microscopic model based on 
relativistic kinetic theory gives a more detailed physical description. Kinetic 
theory ofi^ers a microscopic approach to describe the macroscopic features 
of matter, rather than the phenomenological fluid dynamics and its asso- 
ciated thermodynamics. Starting from a microscopic approach, to obtain 
an effective macroscopic description, is the most fascinating feature of this 
theory. The theory is based on a simple function which is called the distri- 
bution function which is a solution of Boltzmann's or Liouville's equation, 
and describes the dynamics of the system at the microscopic level. At the 
macroscopic level, the mass density, flow density, pressure, and the other 
macroscopic quantities arc obtained from the distribution function. On the 
other hand, observations at large scales, such as stellar systems, can help 
us to improve our understanding of microphysics. 

In chapter 3, general relativistic Liouville's equation in the post-Newtonian 
approximation was studied. In the static case, the equilibrium, two static 
solutions of the Liouville's equation in this approximation are obtained. 
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These integrals are generalizations of the classical energy, E = -^v + + 
(20^ + iIj)/c^, and angular momentum, li = EijkX^v^ exp{—(f)/c'^). In this 
spirit, the polytropic model, a simple model for a neutron star, was studied. 
Our results show that the post-Newtonian corrections tend to reduce the 
radius of any polytrope. This is a consequence of the fact that the post- 
Newtonian correction is more significant for systems with larger density 

i- 

Linear perturbations of phase space distribution functions was investi- 
gated in chapter 4. We introduced the linearized Liouville-Einstein equa- 
tion in this approximation. We showed that, if the underlying potentials 
are spherically symmetric, the evolution equation is 0(3) symmetric, ie. the 
linearized Liouville-Einstein operator commutes with the angular momen- 
tum operator in phase space, Sijk^x^ + -^). Then the modes can be 
characterized by a pair of angular momentum eigennumbers, {j,m). The 
eigenvalues Uj are, however, (2j + 1) fold degenerate. 

Furthermore, we showed that the post-Newtonian gravitational poten- 
tials may excite some of the neutral modes of the star and that these modes 
are purely relativistic effects. Using the 0(3) property of pnl, we proposed 
distribution functions for perturbations that are functions of classical energy 
and classical angular momentum, Eqs. ( [4.18| ) and ( |D.8D : 



= Re A,^ ; = a/(e, Z^) Ji+"/i = bf{e, 



+ 5 



Although these functions are neutral in classical approximation, they are 
not so in pn order. Neutral, here, means to belong to zero frequency modes. 
The weak pn forces generate a sequence of low frequency modes from such 
perturbations. In their hydrodynamic behavior, they constitute a sequence 
of low frequency toroidal modes. There is an oscillatory goi component of 
the metric tensor associated with these modes. From a conceptual point of 
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view, they are similar to toroidal modes of slowly rotating fluids generated 
by Coriolis forces or to the standing Alfven waves of a weakly magnetized 
fluids §. 

The latter perturbations are analogue of the recent quasi normal modes 
in relativistic systems believed to have been originated from the perturba- 
tions of the space-time metric, gravitational wave modes (ly-modes). Kokko- 
tas and Schutz ||113|| first recognized the w-modes in a toy model of a finite 
string (to mimic a fluid) coupled to a semi-infinite one (to substitute the 
dynamical space-time). Such a system accommodates a family of damped 
oscillations due to the emission of gravitational wave. Different investiga- 
tors have proposed different mathematical and numerical schemes to isolate 
these modes ||114|| - [ll20| . They verified that strongly damped (w-) modes, 
due to the space-time metric perturbation, do indeed exist in realistic stellar 
models. 

In chapters 5 and 6, we studied the recent and interesting instability 
in rotating neutron stars. Recently it has been shown that the instability 
of perturbations of rotating stars are important during the early history 
of hot neutron stars. These perturbations are driven by the Coriolis force 
which is always present in a rotating star and are known as r-modes. The 
instability of these modes cause a rotating neutron star's rotation rate to 
slow down, emitting gravitational radiation in the process. This emission 
of gravitational radiation is important both as a possibly detectable source 
and as a mechanism to explain the observed spin rate of neutron stars. 

It is important to improve our understanding of the various factors that 
go into r-mode stability analysis. One important aspect which has needed 
further elaboration is the role of dissipation in the fluid. Normally, the 
fluid's viscosity is thought to damp out any instability if the star is relatively 
cool. However, the earliest analyses only considered a model for the viscosity 
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based on the Navier-Stokes theory (which is known to have serious problems, 
such as faster-than-hght propagation of signals). 

As a first attempt to include more realistic microphysics, we have consid- 
ered the role of vorticity on the stability of r-modes. This effect is predicted 
by kinetic theory when the unperturbed equilibrium state is rotating, but 
is absent in Navier-Stokes theory. In standard Navier-Stokes theory, the 
angular velocity of the fluid has no effect on viscous stress or heat flux. We 
calculated the vorticity- shear viscosity coupling and showed that the cou- 
pling between vorticity and shear viscous stress predicted by kinetic theory 
can in principle have a significant effect on r-mode instability in neutron 
stars. The Miiller-Israel-Stewart correction of Navier-Stokes theory predicts 
that colder stars can remain stable at higher spin rates, so that accreting 
spin-up could be protected from r-mode instability 0]. 

Normally all neutron stars which have been observed are seen to be 
rotating while showing a slow down of their rotation rate. However, some 
pulsars, such as the Crab pulsar exhibit glitches, which are brief periods 
during which their rotation rate suddenly increases. We have studied the 
role of r-modes in the post-glitch relaxation of radio pulsars. We have 
shown that excitation of the r-modes at a glitch may provide a solution to 
an unsolved observed effect in post-glitch relaxation of the Crab pulsar [Q. 

Of course, our analysis is limited by the fact that we have followed the 
standard assumption in viscous stability analysis and ignored superfiuid ef- 
fects that will become important at lower temperatures (see, e.g., ||111|| ). 
Superfiuid "friction" effects are thought to prevent /-mode instability, and 
these effects are likely to be relevant also for r-modes. These effects may 
strongly alter the vorticity correction effect, and the possibility of r-modes 
excitation in the Crab like pulsars (T ~ 10® K). 
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Appendix A 

The post-Newtonian 
approximation 



The Einstein field equations are nonlinear, and therefore cannot in general 
be solved exactly. In most cases, by imposing some symmetries such as time 
independence, spatial isotropy and/or homogenity, we were able to find some 
exact solutions, the Schwarzschild and the Freidmann-Robertson- Walker 
metrics for examples. But we cannot actually make use of the symmetries 
in all problems. Solar system is the familiar example of non-static and 
anisotropic case. 

In most problems what we need is not to find the exact solutions of the 
problems, but we need a systematic approximation method to extract the 
solutions without any assumed symmetry properties of the problem. 



The post-Newtonian approximation was historically derived [[132| , |133 



134|| , to the study of the problem of motion. But in the last three decades. 
It is used largely to study dynamics of stellar systems like compact stars 
and black holes. 

In this appendix we introduce the post-Newtonain approximation in 
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details. We follow |^2| to present this method. 

Consider a system of particles that, like the sun and the planets, are 
bound together by their mutual attraction. Let M, f, and v be typical 
values of the masses, separations, and velocities of these particles. The 
Newtonian typical kinetic energy Mv^ j1 will be roughly of the same order 
of magnitude as the typical potential energy GM jf, so 

2 GM 

^ ^— . A.l 

r 

The relation will be exact for a particle moves with velocity v in circular or- 
bit of radius f about a central mass M. The post-Newtonian approximation 
may be described as a method for obtaining the motion of the system to 
one higher power of the small parameters GM/f and v'^ than given by New- 
tonian mechanics. It is also referred to as an expansion in inverse powers of 
the speed of light, c. Here we prefer to use v/c as the expansion parameter. 

From our experience with the Schwarzschild solution, we expect that it 
should be possible to find a coordinate system in which the metric tensor is 
nearly equal to the Minkowski tensor r/^j^, the corrections being expandable 
in powers of v/c. In particular, we expect 

goo = -1 + ^5'oo + ^goo H , (A. 2a) 

gij = Sij + '^gij + ^Qij H , (A.2b) 

goi = ^goi + ^goi + ■■■, (A.2c) 

where the symbol ^gfj,u denotes the term in g^i, of order (v/c)^ . Odd powers 
of v/c occur in gio because gio must change sign under time- reversal trans- 
formation t —t. These expansion lead to a consistent solution of Einstein 
field equations. The inverse of the metric tensor is defined by the equation 

g'^go^. = g'°goo + g''goj = o, (A.3a) 
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g'^'go^. - g'^'goo + g'^goi = i, (A.3b) 

g'^'gj,. = g'^gjo + g'^gjk = Sij. (A.Sc) 

We expect that 

g''^5,,+ V+ V + (A.4b) 

^oi = + + • • • , (A.4c) 

and inserting these expansions into the Eqs. (A. 3), we find 

= - '^oo; ^g'' = - ^g'^ = 'goi- (A.5) 

The affine connection may be obtained from the famihar formula 



In computing F^^ we must note that the scales of distance and time, f and 
fjv^ respectively. So the space and time derivatives should be regarded as 
being of order 

d \ d vie 
dx^ f' cdt f 
Using the metric expansions we find that the various components of F^^^ 
have the expansions 

ri. = + + • • • ; (fo^ r^o, rjk, f°,), (A.7a) 
ri. = 'r;:, + 'rj, + • • • ; (for ri-, f°„, F°j). (A.7b) 

The symbol ^F^^^ denoting the term in F^^^ of order (v/c)^ /f. After some 
manipulating one finds 

=rL=-i^^. (A.8a) 
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4T-ii 



oj 



2pi 



3t-io 



2t-io 



2 (9a;* 

3, 



cat 

2, 



+ 2 ^^-'"ax^- 



1 / d 9oi d ^Qoi d^gjo ^ 

2 V dx^ cdt dx^ y 

2 \ dx'' dx^ 



dx^ 

1 (d'go. 

2 V cdt 

2 1 



r". = 0. 



The Ricci tensor is defined by 



p — pA 



-pA _ -nri pA 



(A. 8b) 
(A.8c) 

(A.Sd) 

(A.Se) 

(A.Sf) 
(A.Sg) 

(A.9) 



Using Eqs. (A.6)-(A.7) we find that the components of R^i, have the ex- 
pansions 



Rqo — -^00 ~t~ -^00 ~t~ 



Inserting Eqs. (A. 6) in Eqs. (A.9) we obtain 

Q2pi 



Rn 



dx^ 



d ^n, d 



cdt 



3p _ 

cdt 



dx^ 
aspi. 

dx^ 



_l_ 2po 2pi 2pt 2-pj 

oi 00 oo ijy 



'R 



tj 



dx^ ' dx^ dx^ 
Therefore in terms of metric tensor, the Ricci tensor will be 
1. 



'R 



00 - ^9o 



(A.lOa) 
(A.lOb) 
(A. 10c) 

(A.lla) 
(A.llb) 
(A.llc) 
(A.lld) 

(A.12a) 
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1Q2 2 a2 3 1 1 q2 2 



00 



1 /a / d '^9oo\ , ^ (d ^goo\ ( d 



+ 



2 \ dx^ J \ dx"^ J A\dx^J\ dx^ J 

-\ W) (^) ' 

" 2 c^x^-^t 2 dx^dxJ 2 c^x^St ^ 2 ^ ^ ^ 

20 _ I g'Voo , I g^Vfc I g^V, 1^22^ /A-,o,x 

~~ 2 ^ 2 2 ^x'^^x^ ^ 2 ^'^^ ^ ^ ^ 

By choosing a suitable coordinates system, one can simplify the above equa- 
tions. It is always possible to define the x^ so that they obey the harmonic 
conditions 

r^l. = 0. (A.13) 

Substituting Eqs. (A. 3) and (A. 7) in Eq. (A. 12), we find that the vanishing 
of the third-order term in g^'^T"^^ gives 

1 d ^goo d ^goi 1 d ^gu / a 1 q ^ 

while the vanishing of the second order term in gives 

ld^g^^dp^_ld^^^ (A.lSb) 



2 dx' dx^ 2 dx' 



But 



cdt ^''^ 2 c^dt^ cdx^dt 2 c^dt^ 

dx^ ^^-J cdt cdtdx^ dx^dxi 

-1^=0 
2 cdtdx^ 



Q^k \^ f^'^j Q^i f^'-j dxWx^ dx'dxi dx^dx^ 

2o 

dx'^dx^ 
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So Eqs. (A. 11) now give simplified formulas for the Ricci tensor 

'Roo = '9oo, (A.14a) 

4 D V72 4 yoo -'-2 yoo , J- 2 \2 / a i ^i \ 

Roo = 2 V - - - g.,^-^ + -(V , (A.14b) 

'i?oi = ^Vyoi, (A. 14c) 

'R,j = (A.14d) 

The Einstein field equations are 

R,^ = -^-^{T,u-l9,.T^). (A.15) 

The various components of energy momentum tensor will have the expan- 
sions 

rpOO ^ OrpOO ^ 2rpOO ^ . . . ^ (A.16a) 

j^ij _ 2j.ij ^ 4j.ij ^ . . . ^ (A. 16b) 

r°^= ^T°'+ 3t°^ + --., (A. 16c) 

where -^T''^ denotes the term in T^^^ of order {M /f^){v/c)^ . Therefore 

— '^ni' ~ 29iJ.vTx, (A. 17) 

we find 

500 — " Soo + Sqo + ■ ■ ■ ) (A. 17a) 
Sij^ "Si^+ + (A.17b) 

501 — + ^Soi + • • • . (A. 17c) 



In particular 



''Soo = ^ "T'"', (A.lSa) 
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^Soo = - ( ^T°° - 2 ^goo "T°° + ^T') , (A. 18b) 

"Sij = \ °T''%j, (A.18c) 

i^^i = - iT''^ (A.lSd) 

Using Eqs. (A. 14) and (A. 18) in field equations, we find that the field 
equations in harmonic coordinates conditions 

'goo = (A.19a) 

SttG ^ 2y„„ _ 2 2^^^ ^ 2ynj ^ (A.19b) 



c4 

'go^ = ^ ^T-, (A.19C) 

2^,, = "T-. (A.19d) 

Prom Eq. (A. 19a) we find 

""goo = -20, (A.20a) 
where is the Newtonian potential, defined by Poisson's equation 



= ^ (A.20b) 



Also ^g'"" must vanish at infinity, so the solution is 

G f ''T°°(x',t) 



0(x,t) = -^/j^^^rfV. (A.20C) 

Prom Eq. (A.19d) we find that the solution for "^g^j that vanishes at infinity 
is 

'gij = -205,,-. (A.21) 

^gio is a new vector potential 

^goi = ^i, (A.22a) 
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and the solution of Eq. (A. 19c) that vanishes at infinity is 

6(x,i) = -AG J I 1^ "^'- (A-22b) 

Using Eqs. (A. 19b) and (A.20) and the identity 

we obtain 

V = -20' - 2^. (A.23a) 
The scalar potential t/j satisfies 

= ^ + 471^ f ^T"" +2 r*) , (A.23b) 
at ^ ' 

with solution 

^(x, t) = -4 / (i^^^ + G ^^"°(-'' ^) + G ^T-i^, t)) . 

(A.23c) 

The coordinates condition, Eqs. (A. 13), imposes on (p and ^ the further 
relation 

f + V.f = 0. (A.24) 

while the other coordinate condition, Eq. (A. 12b) , is automatically satisfied. 
The various components of the affine connection are 

^ (A.25a) 
U=|,(20V^) + §. (A.26b) 

»C ^ ^, (A.25e) 
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^|.V*. (A.25i) 
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Appendix B 
Derivation of Eqs. (3.5) 



Consider a general coordinate transformation {X, U) — {X^, to (Y, V) — 
{Y^, V^). The corresponding partial derivatives transform as 



d/dX 
d/dU 



M 
/ 

v 



d/dY 
d/dV 

dY/dX dV/dX 

dY/dU dV/dU 



' d/dY 
^ d/dV 



(B.l) 

where M is the 7x7 Jacobian matrix of transformation. Setting X — Y — 

x^, V = v"^ and U = W for our problem, one finds 

M^f"-^'^"-^ "'^^"■^ 1, (B.2a) 



dx^'/dx" dv'/dx'' 
dxi'/dW dv'/dU^ 



and 



dx^ldx"" dU^/dx" 
dx^jdv^ dU'/dv^ 



One easily finds 

dxi'/dx'' = 5^^; dxf'/dv^ = 0, 



(B.2b) 



(B.3a) 
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dU'/dx'' = v'dU'^/dx'' = ^° ""'^^f (B.3b) 
duydv' = U°5ij + v'duydv' = U'^Sij - U'^^'gjpv^ (B.3c) 

Inserting the latter in and inverting the result one arrives at M from 
which Eqs. (3.5) can be read out. ae 
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Appendix C 

Post-Newtonian 
hydro dynamics 



Mathematical manipulations in composing of this work has been tasking. To 
ensure that no error has crept in the course of calculations we try to derive 
the post-Newtonian hydrodynamical equations from the post-Newtonian 
Liouville equation derived earlier. From Eq. (3.6a) one has 

= + + ^v2)£^' + a'^]F, (C.l) 

where CJ^ and CF^ are given by Eq. (3.10). We integrate CF^uF over the 
v-space: 

y CrjjFd?v = J [(c' + + ^v')£'^' + JCr]Fd\. (C.2) 
Using Eqs. (3.12) and (3.13), one finds the continuity equation 

^ / Ot-iOO I 2rpQ0\ I / IrpOj I 3rp0j\ _ QrpOO ^4' _ p, 

Cdt^ ^ '^dx^^ ^ + ^ J ^Q^ 

(C.3) 

which is the pn expansion of the continuity equation 

T% = 0, (C.4) 
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Next, we multiply ^^"[/F by f * and integrate over the v-space: 



(C.5) 



After some calculations one finds 

f IrpOi I 3rp0i\ I 9 / 2rpij . 4rp 



cdt 

_|_ OjlOO 



d 



dx^ 



45,; 



90\ 

c9t J 



cdt I 



c^dx^ 



' dx 



dx^ 



(C.6) 



The latter, the pn expansion of 



=0; z = l,2,3. 



(C.7) 



is the same as that of Weinberg pf, QED 
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Appendix D 

Eigensolutions of and Jz 



As pointed out earlier, Jj's of Eq. (4.15) have the angular momentum 
algebra, 

[Jj, Jj] = iSijkJk- (D-1) 

Therefore, the simultaneous eigensolutions of and Jz, Ajj„(x, u), obey 
the following 

J^Ajm = j{j + l)A,m, J = 0, 1, • • • , (D.2) 

JzAjra = mAjm, -j <m< j. (D.3) 
The ladder operators, J± — Jx^iJy, raise and lower the m values: 

J±Ajv„ = ^{]Tm){]±m + l)Kjm±i- (D.4) 

In particular 

J±^j,±j = 0. (D.4a) 
The effect of Ji on classical energy integral, e — /2 — 9{r), and the classical 
angular momentum integral, /j = EijkXjUk, are as follows 

J,e^Jil^^Jif{e,l^)^0, (D.5a) 
Jilj = iSijklk- (D.5b) 
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Theorem 1: 

A,,^, = ii = {\y{i,±iiyy. (D.6) 

Proof: 

JJi = = ±jli, by (D.5b), (D.7a) 

J% = (J_ J+ + + J,)li = j{j + by (D.4a) and (D.7a)j(D.7b) 
jHi = (J+J_ + J,^ - J.)F: = j{j + (D.7c) 

QED. Combining Eqs. (D.6), (D.4) and (D.5) one obtains 

Ajm = af{e, f)Ji^^V_ = 6/(e, f)J^^l^, (D.8) 

where /(e, is an arbitrary function of its arguments, and a and h are nor- 
malization constants. Examples: Aside from an arbitrary factor of classical 
constants of motion, one has 



Aio = ^z, (D.9a) 

Ai±i = i±, (D.9b) 

A2o = 2M_-/^ = ^(3/^-/2), (D.9c) 

A2±i = /±/., (D.9d) 

A2±2 = ^|. (D.9e) 



Theorem 2: The vector field V-'™' = / Aj„jU(ir2 is a toroidal vector field 
belonging to the spherical harmonic numbers (j, m) , where integration is 
over the directions of u. 

Preliminaries: Let (■»?, (/?) and (a, (5) denote the polar angles of x, of u, 
respectively, and 7 be the angle between (x, u). Also choose magnitudes of 
X and u to be unity, for only integrations over the direction angles are of 



120 



concern. One has cos 7 = cos'j? cos a + sini? sinct cos (99 — (3) 

Ur — cos 7, (D.lOa) 

u.ff = — sini? cos a + cos'j? sin a cos((/3 — /?), (D.lOb) 

Uy, = — sin a sm{ip — , (D.lOc) 

1+ = i{sm^ cos a e''^ - cosi? sina e^^). (D.lOd) 

Proof: By induction, we show that a) V-'-' is a toroidal field and b) if V-^"* 
is a toroidal field, so is V^ "^~^. 

a) Direct integrations over a and /3 gives 

Vj^ = j l{urdn - 0, dn^ sin a da dp, (D.lla) 
Vr = J liu,dn = -^-^Yni^, (D-llb) 
V^^ = I liu^dn = g;^yjM V)- QED. (D.llc) 

b) Suppose V-'"* is a toroidal vector field and calculate V-' = / {J^Ajm)ndil, 
where J± = L± + K±, L± = ie^^^^ ± ^cotg^^), K± = ±e^'^{£ ± 
icotga^). Again direct integrations gives 

yj m-1 ^ j^_yjm ^ yjm ^ (D.12a) 

Vr-'^—^.^Yj^.,i^,^), ifVi'^^—^J-YjM^), (D.12b) 

sm V dip sm v dip 



Vr~'--^Yj^_,{^,ip), if V^"^^-^Y,^{^,ip). (D.12C) 

QED. 
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